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EXECUTIVE SUMMARY 


The goal of this project was to develop and test a coupled model-physics suite to 
investigate the interactions between the surface fluxes, the diurnal cycle, the planetary boundary 
layer characteristics and the moisture cycle (including shallow and deep clouds) in a 
nonhydrostatic three-dimensional (3-D) mesoscale model (the Penn State University / National 
Center for Atmospheric Research (PSU/NCAR) mesoscale model, widely know as MM5). 

Atmospheric modeling over continents is incomplete without adequate treatment of the 
land-surface boundary, the planetary boundary layer and its associated cloud systems. It was 
hypothesized that an improved treatment of the regional water cycle can be achieved by using a 
3-D mesoscale numerical model having high-quality parameterizations for the key physical 
processes controlling the water cycle. These included a detailed land-surface parameterization 
(the Parameterization for Land-Atmosphere-Cloud Exchange [PLACE] sub-model of Wetzel and 
Boone, 1995), an advanced boundary- layer parameterization (the 1. 5-order turbulent kinetic 
energy [TKE] predictive scheme of Shafran et al., 2000), and a more complete shallow 
convection parameterization (the hybrid-closure scheme of Deng et ah, 2002a,b) than are 
available in most current models. PLACE is a product of researchers working at NASA's 
Goddard Space Flight Center in Greenbelt, MD. The TKE, shallow convection and deep 
convection schemes are the result of model development at Penn State. 

The long-range goal is to develop an integrated suite of physical sub-models that can be 
used for regional and perhaps global climate studies of the water budget. Therefore, the work 
plan focused on integrating, improving, and testing these parameterizations in the MM5 and 
applying them to study these processes over the Southern Great Plains (SGP). These schemes 
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have been tested extensively through the course of this study and the latter two have been 
improved significantly as a consequence. 

The work conducted in this study provided the first insertion of the PLACE scheme into 
the MM5v3 model (the latest version of the MM5). The result is a fully coupled comprehensive 
parameterization of land-surface processes (PLACE) in the MM5 mesoscale model. Care has 
been taken to ensure appropriate initialization of the surface fields and definition of parameters. 
An offline PLACE system forced by observations was constructed and refined for this study to 
initialize 3-D MM5-PLACE. Extensive analysis of model results was performed against ARM- 
CART and Oklahoma Mesonet data. Utilization of MM5 with PLACE for a summer convective 
case has shown that the coupled model better captured the spatial distribution of moisture in the 
boundary layer. The observed west-east gradient in afternoon surface-layer moisture across 
Oklahoma was reproduced by the PLACE-MM5 runs, but not by the simpler force-restore Slab 
model within MM5. Based on fine-scale data for the subgrid soil-moisture heterogeneity, we 
have shown that it has a positive impact in the computation of surface fluxes of water vapor and 
sensible heat. 

It was found through extensive testing in 3-D and 1-D frameworks that the PLACE land- 
surface scheme performed at about the same level of accuracy as other land-surface schemes 
used to represent the surface fluxes over the Great Plains in normal operations. In some cases, 
PLACE was able to perform better than alternative schemes. However, in many situations, it 
was found that a more simple force-restore land-surface contained originally in the Blackadar 
boundary layer parameterization (Zhang and Anthes 1982), known as the Blackadar Slab 
scheme, generally provided about the same level of skill or better. Thus, while the results of 
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MM5-PLACE simulations were as good or better than any others known to be done with 
PLACE, it proved difficult to show an overall clear advantage when PLACE was used. 

The underlying reason for this result is attributed mostly to the uncertainty involved in 
defining many of the key soil and vegetation parameters that must be defined in order to run 
PLACE. Even with the best soil-characteristic data sets currently available, such as the State 
Soil Geographical Database (STATSGO), PLACE produced soil-moisture evolutions that did not 
always agree well with soil-moisture measurements over the Great Plains. Work is currently 
underway at Penn State to evaluate MM5-PLACE during July 1997 using ESTAR soil moisture 
and NDVI-based leaf-area index and vegetation fraction. 

The 1.5-order TKE-predicting turbulence scheme of Shafran et al. (2000) worked 
reasonably well in dry (unsaturated) conditions, but was found to have deficiencies in saturated 
layers. A reformulation was introduced that now accounts for the effects of saturation on the 
calculations in the buoyancy-production term of the turbulent kinetic energy equation, the 
stability parameters and turbulent length scales. This saturation dependency is missing in many 
current numerical models. The change allows the turbulence sub-model to generate more intense 
mixing in clouds, especially near cloud tops where the longwave radiative flux divergence can 
cause intense cooling. The consequence of this improvement is a more realistic thermodynamic 
and moisture structure in cloudy layers and better representation of fog layers. 

Extensive evaluations and testing of the shallow-cloud parameterization also were 
conducted in the course of this investigation. It was found that this convection sub-model 
contained several layer-dependent formulations for some key quantities related to the cloud-base 
mass flux. In particular, the discrete thickness of the model layers tended to cause a 
discontinuous behavior of the mass-flux calculations that sometimes led to erratic behavior of the 
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convection scheme. This behavior occurred even for vertical configurations that are considered 
to be very good for 3-D models (e.g., layer thicknesses of only -40 m in the boundary layer and 
the cloud layer). Further work led to changes in the formulations of the shallow convection 
scheme so that the dependency on discrete layer thickness was reduced or eliminated. The 
consequence was that the key quantities of the shallow-convection scheme, including mass flux 
and detrainment rates, now evolve in a generally smooth manner as the clouds develop and grow. 
Not only does the reformulation provide a more realistic shallow-cloud structure, but it also 
makes the evolution of calculated water-budget terms smoother and much easier to interpret. 

Following the extensive testing of each individual parameterization, several combined 
experiments were conducted in the 3-D MM5 using PLACE or the Blackadar Slab land-surface 
scheme, and the improved versions of the PSU 1.5-order TKE-predicting turbulence scheme and 
shallow-convection scheme. These experiments showed that the MM5 can produce tairly 
realistic distributions of shallow convective cloud fields in stratocumulus and fair-weather 
cumulus (humilis) conditions. Cloud area and liquid water content are reasonable and are within 
the accuracy expected, given typical uncertainties in the initial atmospheric and soil states. The 
vertical structure for the resolved-scale temperature and moisture fields in the cloud layer appear 
to evolve in a generally realistic way, although further testing is required. 

The coupling of the PLACE-MM5 model with the new PSU shallow convection scheme 
in MM5 and applied to the same case has shown that we are able to include in our simulations 
important physical interactions that control the development of convection over land (cloud- 
radiative interactions, cloud-surface interactions, and cloud-boundary layer interactions). This 
evaluation of the complete 3-D modeling system in a complex convective situation over land has 
also suggested several improvements that are being tested in the new parameterizations, 
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especially in the shallow convection and turbulence sub-models. We have shown that, as we are 
able to incorporate more physically realistic interactions in our model system, the accuracy of 
each sub-model must be increased. Otherwise there is a possibility that the simulations may 
diverge from the observations even more rapidly than without the inclusion of these physical 
processes. However, our results also suggest that the inclusion of these complex physical 
interactions is crucial if one wants to better simulate the development and evolution of 
convection over land. 
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1 . 


INTRODUCTION 


1.1 Objective 

The objective of this investigation was to study the role of shallow convection on the 
regional water cycle of the Mississippi and Little Washita Basins using a 3-D mesoscale model, 
the PSU/NCAR MM5. The underlying premise of the project was that current modeling of 
regional-scale climate and moisture cycles over the continents is deficient without adequate 
treatment of shallow convection. 

It was hypothesized that an improved treatment of the regional water cycle can be 
achieved by using a 3-D mesoscale numerical model having a detailed land-surface 
parameterization, an advanced boundary-layer parameterization, and a more complete shallow 
convection parameterization than are available in most current models. The methodology was 
based on the application in the MM5 of new or recently improved parameterizations covering 
these three physical processes. Therefore, the work plan focused on integrating, improving, and 
testing these parameterizations in the MM5 and applying them to study water-cycle processes 
over the Southern Great Plains (SGP): 

(1) the Parameterization for Land- Atmosphere-Cloud Exchange (PLACE) described by 
Wetzel and Boone (1995), 

(2) the 1.5 -order turbulent kinetic energy [TKE]-predicting scheme of Shafran et al. 
(2000), and 

(3) the hybrid-closure sub-grid shallow convection parameterization of Deng (1999). 
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Each of these schemes has been tested extensively through this study and the latter two have 
been improved significantly to extend and their capabilities. 

1.2 Background 

The summertime regional-scale water budget over continents is an important factor in the 
overall global water cycle. Although horizontal transport in numerical models is relatively 
straightforward, improvement of vertical transports, especially for moisture, represents a crucial 
need. Summer regimes in middle latitudes generally do not exhibit the strong dynamical forcing 
and large-scale vertical motion fields common in winter cyclones. Instead, vertical moisture 
transport is highly dependent on radiation-driven sub-grid scale physical processes (Dickinson 
1995). These include land-surface fluxes, boundary layer turbulent processes, and both shallow 
and deep convection. 

Modeling of regional-scale climate and moisture cycles is incomplete without adequate 
treatment of shallow convection (Dickinson 1995, Wetzel and Boone 1995), which strongly 
affects the surface fluxes, the upper boundary of the mixed layer and the development of 
convective precipitation. Shallow convection is defined here as the process by which air parcels 
become saturated, such that non-precipitating, buoyancy-driven clouds form, dissipate, and 
interact with their environment. Although shallow convective clouds can form at any level in the 
atmosphere, our attention is confined to the layer near the top of the planetary boundary layer, 
where these clouds are most often observed. Shallow convection includes both cumulus and 
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stratocumulus clouds types of various depths and area coverage, and may involve both water and 
ice phases. 

Despite the acknowledged importance of shallow convection in moisture transport and 
climate studies, it is often treated in a relatively simply manner in regional and global-scale 
numerical models. For example, Betts and Miller (1986) employ a simple adjustment scheme 
that drives the environment toward a pre-determined state, but does not provide many details of 
the cloud structure. In addition, their approach for shallow convection appears to be best adapted 
to the tropical marine environment. Seth and Giorgi (1996) studied how organized mesoscale 
circulations induced by vegetation can affect vertical transports of heat and moisture. Using a 
regional-scale climate-adapted version of the hydrostatic MM4 model (RegCM2), they found 
important effects on monthly precipitation resulting from land-surface heterogeneity. However, 
their model treats the effects of shallow clouds only in a very crude sense. This could be a 
serious omission because of the strong influence these clouds have on the radiation budget and 
vertical moisture transport. 

Although shallow convection produces very little (if any) precipitation and does not 
affect its environment nearly as much or as rapidly as deep precipitating convection, its effects 
are not trivial. Shallow convection can be important for altering local and regional-scale 
precipitation patterns and the structure of temperature, water vapor, and mixing depth in the 
turbulent boundary layer. In a similar way, shallow-cloud patterns and their area coverage are 
affected by horizontal heterogeneity in the land-surface processes, even over fairly small areas 
(Wetzel and Boone 1995, Wetzel et al. 1996). Shallow convection at the top of the mixed layer 
can act as a sink for boundary-layer moisture, while transporting water, aerosols and chemical 
species into the middle troposphere. Shallow clouds also have an impact on subsequent 
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convective precipitation by changing the environmental temperature and moisture structure in 
which deep clouds must grow. By using up some of the convective available potential energy 
(CAPE) of boundary-layer air, shallow convection can delay or in some cases even prevent the 
onset of deep convection. In other cases convection could occur sooner due to pre-conditioning 
of the environment. 

1.3 Organization 

Section 2 documents the shallow-convection parameterization, which is the central sub- 
model needed for this investigation. This parameterization, developed by Deng (1999) is 
currently under review for publication (Deng et al. 2002a,b), so its design is presented here to 
assist in understanding the work carried out through this study. The PLACE land-surface 
scheme and the TKE-predicting turbulence scheme already have been described in the literature 
by Wetzel and Boone (1995) and Shafran et al. (2000), so it is unnecessary to repeat those 
descriptions here. 

Next, Section 3 describes the most important 1-D experimentation conducted in this 
study, using a single column model (SCM) based on the vertical structure of the MM5 model, 
and the offline (stand-alone) PLACE methodology used to initialize the 3-D MM5 PLACE. This 
1-D framework proved very helpful for understanding and refining the parameterizations for 
turbulence and shallow convection. Section 4 describes the 3-D numerical experiments carried 
out in the MM5 and discusses the implications of the model results. Finally, Section 5 
summarizes what has been learned through this study and some of the most important work 
remains for the future. 
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2 . 


THE SHALLOW CONVECTION PARAMETERIZATION 


The Kain-Fritsch (1990) deep-convection scheme is the progenitor for the present 
shallow-cloud development. In the Kain-Fritsch (KF) scheme, if a cloud updraft fails to reach a 
critical depth necessary to support rain (i.e., shallow cloud only), it is assumed to disappear 
without having any impact at all on its environment. Here, we seek especially to describe how 
those shallow clouds grow and interact with their surroundings. The major elements of the 
shallow-cloud scheme are: 


• definition of initial cloud-parcel characteristics and the convective trigger 
mechanism; 

• a convective-cloud sub-model based on the parcel buoyancy equation, 

• closure assumptions that determine the cloud-base mass flux; and 

• a prognostic scheme for the area and water content of clouds that result from 
the detrainment of the convective-updraft air. 

The cloud parameterization is designed to represent the physical linkage between the 
turbulent planetary boundary layer (PBL) and moist convective processes in a multi-layer 
framework. We define the PBL as the layer beginning at the surface and extending upward to 
the point where surface-based turbulence rapidly decreases with height. Although the top of this 
layer can be identified in a model by the drop-off in turbulent kinetic energy, it also tends to be 
characterized by a sharp jump in the potential temperature. Thus, it can be distinguished from 
the cloud-topped boundary layer (CTBL) discussed by some investigators (e.g., Albrecht 1979, 
Agee 1987), which includes the PBL plus a shallow-cloud layer having distinctly different 

stability characteristics. 

We assume that shallow convective clouds consist of active updrafts and (approximately) 
neutrally buoyant clouds (NBCs). The NBCs represent either the remnants of previous updrafts 
or the cloud mass dispelled from currently active updraft cores. This production of NBCs from 
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detrained updraft air is consistent with the conceptual model of shallow convection presented by 
Wyant et al. (1997). The most vigorous convective updrafts generally cover less than 10% of a 
grid area, while the associated NBCs may cover a much greater fraction of the sky. Although the 
shallow updrafts are assumed to produce no rain, the NBCs may generate some light 
precipitation. In general, most of the updraft mass is detrained near the cloud tops as it reaches 
or overshoots the equilibrium level. Over time, this detrainment process and the consequent 
induced subsidence can change or even dominate the thermal and moisture structure of the 
mesoscale environment. If the shallow-cloud updraft exceeds a critical depth (defined as Dkf - 
4 km in this application), it is considered to transition to deep convection (thunderstorm), often 
with heavy precipitation and strong moist downdrafts (Kain and Fritsch 1990). On the other 
hand, under conditions with a strong capping inversion and large vertical moisture flux, the 
detrainment process can lead to an accumulation of vapor and detrained cloud mass at the 
inversion base, so that a solid stratus deck may develop. Thus, the shallow convection can have 
a direct link to both resolved stratiform cloud and sub-grid deep convective cloud, both of which 
can be active precipitation generators. In many 3-D models, representation of the physical 
relationships among these three types of cloud is weak at best. 

2.1 Cloud-Parcel Initial Characteristics 

The parameterization is built as a one-dimensional column sub-model. Active cloud 
updrafts are triggered when parcels originating in the PBL are able to reach their lifting 
condensation level (LCL). The characteristics of a potential cloud-initiating parcel are defined 
by its virtual potential temperature, 6 vp , water-vapor mixing ratio, q vp , and vertical velocity. The 
values of and q vp for the parcel are defined at each time step from the average ambient values 
of the model layers in the lowest 20% of the PBL or in the lowest two model layers, whichever is 
deeper. The underlying assumptions are that energetic turbulent eddies are initiated near the 
surface and that the largest (and most buoyant) of these rise with only modest dilution through 
the entire PBL to approach or reach their LCL. The LCL for these parcels is calculated 
according to the method of Fritsch and Chappell (1980) and may lie below or above the PBL top. 
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In theory, an ensemble of shallow clouds could be created by defining different parcels with 
initial characteristics based on a probability distribution about 0 V p and (]vpi but that approach has 
not been introduced here. 

Following eddy-transport theory, we hypothesize that these cloud-forming parcels have 
positively correlated thermal, moisture and vertical velocity perturbations. The eddy vertical 
velocity, wy, is scaled to the maximum turbulent kinetic energy (TKE) in the PBL according to 



where the constant ^ results from the TKE definition, assuming that the turbulence is isotropic. 

We then define the release height of the cloud-initiating parcel, zr, to be the lower of two levels: 
the top of the PBL ( h PBL ) or the LCL. The total vertical velocity of the cloud-initiating parcel is 
estimated as Wp — w + wy , where w is the resolved-scale vertical motion at zr • The TKE is 
calculated using a 1.5-order turbulence scheme (Gayno 1994, Shafran et al. 2000). Typically, in 
a weakly forced marine boundary layer w T is only ~0. 1-0.3 ms 1 . Over land in a shear-driven 

nocturnal boundary layer, w T may grow to -0.3-0.8 ms' 1 , while during the afternoon in a 

convectively unstable boundary layer it can easily reach -1-2 ms 1 . If a cloud-initiating parcel 
fails to reach the LCL, it is assumed to return to its point of origin in the PBL without affecting 
the environment. 

Additional factors (e.g., sub-grid terrain irregularity or land-use variability) also may 
contribute to the parcel vertical velocity in certain cases, but most of these are ignored at present. 
However, nonhydrostatic pressure-gradient forces beneath strong deep-cloud updrafts 
(Schlesinger 1984), w N n, are considered when the cloud depth, D c , becomes greater than the 
critical Kain-Fritsch depth, D KF - 4 km. In this case, w NH is defined as a simple function of 
w(z 10 oo) r_1 , which is the parcel vertical velocity at 1 km above the cloud base from the previous 
time step, according to 
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(2.2) 


_ Jo. 15 • [w(z 1000 ) T ' ~Wp] if wUiooo) 7 ' - w p >0and D c > D KF 

W^// — 4 

[ 0 otherwise 

The parameter 0.15 in (2.2) yields a magnitude of w NH in the range 0.5 <w NH < 1.5 ms' 1 when 
deep convection is present. The total parcel vertical velocity at z R , therefore, is given by: 

w R =w + w T +w NH =w p + w NH (2.3) 


If the cloud-initiating parcel reaches the LCL, a convective-cloud updraft forms. Its 
subsequent acceleration and mixing with the environment are calculated using the entraining- 
detraining cloud model of Kain and Fritsch (1990). The entrainment rate at any level is a 
function of the radius of the updraft, R c , and the local parcel buoyancy with respect to its 
environment. While the Kain-Fritsch deep-convection scheme assumes a constant updraft radius 
of 1.5 km, the shallow-convection scheme allows the radius R c to grow smoothly from a 

minimum of R r =0.15 km to a maximum of R Cmax = 1.50 km, which occurs when the 

t min v- l,kUA 

shallow convection transitions to deep precipitating convection. 

For our purposes, we assume that the most significant factors controlling R c are the 

depth of the boundary layer, h PBL , and the depth of the cloud (D c = z T - z B , where z T is the 
height of cloud top and z B - LCL is the cloud-base height). The calculation of the cloud depth 
is described below in Section 2.2. For simplicity all updrafts in a grid cell are considered to have 
equal radius and depth at a given time. This uniform geometry is convenient, although an 
ensemble of different cloud sizes would be more realistic (e.g., Siebesma and Cuijpers 1995). 
The dependency of R c on boundary layer depth is assumed to be greatest for comparatively 

shallow cloud depths and shallow boundary-layer depths. These basic assumptions are 
approximated in the following relationship for the cloud-updraft radius (when D c < D KF ). 


Rr = 




PBL 


and R, 


C min 


< /? c < R c ma X 


(2.4) 
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where b = 1 + 2r J h PBL rj= 12D c_ ; and Ro D and h PBL are in km. Figure 2.1.1 
2 4 -D c 

demonstrates the relationship between R c , D c and h PBL given by (2.4). For deep boundary 
layers notice that even fairly shallow clouds have large radii, while for very shallow boundary 
layers, the clouds must become quite deep before large radii develop. When D c > D KF , 
R c - 1.5km and thermodynamic control is passed to the KF deep convection scheme, which 

includes convective rainfall and moist downdraft development. Obviously, the specific values 
used in this formulation are somewhat arbitrary, but the distributions appear to be qualitatively 
consistent with reality. The application of R c was most effective when averaged over 2-3 time 

steps. 
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2.2 The Convective-Cloud Sub-Model 


As mentioned in Section 2.1, the KF 1-D entraining-detraining cloud model is used to 
calculate the updraft vertical velocity and liquid-water profiles based on the parcel buoyancy 
equation. Unlike KF, however, clouds are not assumed to grow instantly to their mature 
equilibrium level as soon as a cloud is triggered. Here, the cloud top grows at a rate proportional 
to the maximum vertical velocity of the updraft, VF max , estimated to be 

^ = 0.2 •H'„, (2-5) 

dt 

until the equilibrium level is reached. The empirical factor 0.2 imposes a reasonable cloud 
growth rate (e.g., see Simpson 1983), so that under deep-convection conditions, it usually takes 
20-30 minutes to reach the tropopause from the level of free convection (LFC). This approach 
allows detrainment from the growing cloud top to moisten the environment, which is 
oversimplified in the instantaneous cloud growth of the KF scheme. The interpretation of the 0.2 
factor is that cloud-top growth is slowed because the updraft must push aside environmental air 
to continue its upward progress, while the parcel buoyancy equation merely describes the cloud's 
velocity profile in its full-grown state. 

Although the thermal and moisture characteristics of the cloud-forming parcel are defined 
from the lowest 20% of the PBL (Section 2.1), its mass is taken directly from the sub-cloud 
layers nearest to the cloud base. The depth of this updraft source layer, D s , grows as a function 
of the updraft radius from a minimum of 100 m to a maximum of 600 m, according to 

D s =d^l0O0(R c -R Cnan )/d 2 (2.6) 

with the constraint D s < h PBL and where d x = 100 m and dz = 2.7. Thus, as clouds grow wider 

and deeper they are expected to have greater mass flux at cloud base and are likely to entrain air 
from a deeper sub-cloud layer. To satisfy continuity requirements, subsidence is induced in the 
cloud environment to compensate for the mass extracted from the PBL source layer. If the cloud 
remains "shallow", ( D c < D KF ), no sub-grid scale convective downdrafts are allowed to form, so 
all compensation for the upward mass flux must occur through this subsidence mechanism. 
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It should be noted that, although the growth of the updrafts occurs gradually, there is still 
only one size of clouds allowed in a particular grid column in the present formulation. That is, 
unlike the atmosphere, we do not attempt at this time to represent an ensemble of cloud sizes, as 
has been done for deep convection by Arakawa and Schubert (1974). 

2.3 Cloud Updraft Closure Assumptions 

The cloud-base mass flux closure is adopted to determine the intensity of sub-grid 
convection from resolved-scale quantities. The grid-cell mass flux at cloud base is defined as: 

p B =N{nR 2 c pw B ) (2.7) 

where N is the number of updrafts in the cell, w B is the parcel vertical velocity at cloud base, 
p is the parcel air density (here, R c is in m). Generally, closure requires that either p B or N must 
be specified to allow the other quantity to be diagnosed. Of course, for any scheme, errors may 
occur in the calculation of N , R c , or w B , but the key is to estimate p B with reasonable accuracy 
so as to simulate realistic measurable cloud-field characteristics (e.g., cloud fraction, depth, 
liquid-water path length, etc.). Moreover, as stated above, we have made the simplifying 
assumption that all shallow clouds in a grid cell have the same geometry, rather than use a more 
realistic distribution of cloud sizes. Thus, the calculations for N and Rc, in particular, should not 
be considered literal, but merely provide a qualitative estimate needed for the mass flux 
calculations. 

As part of the preliminary development of the updraft module for the shallow-convection 
scheme, five different mass-flux closure assumptions were tested. Brief descriptions and 
comments about these closures are as follows: 

2.3.1 Boundary-Layer Vapor Balance (BLVB) 

Used by Tiedtke (1989) to simulate subtropical trade cumuli, this closure assumes total 
water vapor in the PBL is constant (in the absence of rain). That is, the rate of vapor removal 
from the PBL by cloud-base updrafts is balanced by the sum of surface evaporation and 
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entrainment flux at the PBL top. Thus, A, (the number of cloud updrafts for the BLVB closure) 
can be diagnosed from (2.7) by defining q v )iB as the sum of these two fluxes. Although the 
BLVB closure was found to be reasonable for many marine environments, it greatly under- 
estimated the cloud-base mass flux in continental applications (not shown) and therefore was 
eliminated from further consideration. 

2.3.2 Convective Available Potential Energy (CAPE) Removal 

This closure assumes that total cloud-base mass flux proceeds at a rate necessary to 
stabilize the column over a deep-convective time period, which is generally about 30 minutes 
(Fritsch and Chappell 1980, Kain and Fritsch 1990). Similar to the BLVB closure, the CAPE- 
Removal closure diagnoses the number of updrafts in the grid cell, Aj, from (2.7) using the 
CAPE-derived jUb . The CAPE closure often works poorly in shallow-convection environments 
where the cloud depths are about 1 km or less (i.e., there is little or no CAPE in the shallow- 
cloud layer). Nevertheless, it is still attractive for cases in which most of the shallow cloud lies 
above the LFC. Application of A 2 was found to be most effective when averaged over ~15 time 

steps. 


2.3.3 Boundary-Layer Depth (BED) Relationship 


This closure assumes that the number of updrafts for a shallow-convection environment, 
Nj, is a direct function of the scale of the largest, most-energetic turbulent eddies in the PBL. 
Thus, N 3 depends on the depth of that layer. Since the maximum amplitude of the eddy vertical 
velocity spectra occurs at wavelength ~l.5h PBL (Young 1987), we hypothesize that under 
convective conditions, this geometry can be used to estimate the distribution and number of 
cloud-initiating updrafts in a grid cell according to 


A 


3 


AX AY 
(1.5 h PBL ) 


(2.8) 
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where AXAT is the grid-cell area. The total mass flux, fa, can then be calculated from (2.7). A 
similar closure has been used by other investigators to study boundary-layer rolls (Stull 1988), so 
(2.8) could be modified for use in strongly sheared environments. The BLD closure was tested 
in both continental and marine environments and found to be suitable for fairly shallow cloud 
depths (Dc < 2 km). However, it sometimes caused large oscillations in the mass flux due to 
feedbacks between h PBL , N 3 and fa. 


2.3.4 TKE-Based Closure 


This closure operates on the assumption that shallow clouds basically are driven by the 
TKE in the PBL. Specifically, it scales the shallow-cloud mass flux by the magnitude of the 
maximum diagnosed TKE in the layer. The closure requires that we define four quantities: Ms, 
an upper bound on the amount of mass that can be processed by the shallow clouds (given as the 
total mass in the updraft source layer at the time the scheme is called); t sc , a characteristic 


relaxation time scale; and TKEu, and TKEl, which are upper and lower bounds on the TKE used 
for the mass-flux calculations. The relaxation time scale is defined as T sc = bfD s / w B ), where 

D s is the depth of the updraft source layer, w B , is the cloud-base updraft velocity, and b\ - 4.5 is 
an empirical constant. Then, letting TKEs be the maximum diagnosed TKE in the source layer 
and setting the constraints TKEl, - 1.0 J kg' 1 and TKEu = 10.0 J kg' 1 , we obtain 

TKE l for TKE S < TKE L 

TKE ; = | TKE S for TKE L < TKE S < TKE U 1 , (2.9) 

TKEy for TKE S > TKE U 

where TKEi becomes the effective maximum value of TKE in the source layer. Given these 
reasonable constraints, the maximum mass flux that could occur is given by 

M c 


// fl (max) = 


Tsc 


( 2 . 10 ) 


and the actual mass flux given by the closure becomes 

TKE, 


( 


TKE, 


-)^s(max). 


( 2 . 11 ) 
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For this closure the number of clouds in the grid cell can be recovered from the cloud-base mass 
flux as N 4 = jUb/jUbs, where (2.7) defined the single-cloud mass flux as / i BS = nR c pw B - 

To interpret (2.11), let us describe the "shallow-convection period" to be the time 
required for N 4 updrafts to remove all of the mass from the source layer, M s , (and assuming no 
compensating mass flux into the layer). Then, (2.10) and (2.11) reveal that when TKE, = TKEu 
(i.e., when TKE S > 10.0 J kg' 1 ), the shallow-convection period equals T sc ; however, when TKEi 

< TKEu (the normal condition), mass removal occurs more slowly by the factor TKEi/TKEi. For 
example, in a typical convectively unstable continental environment, we might find that TKE, - 
2.0 J kg' 1 , w B - 1.15 m s' 1 and D s is -300 m. Thus, the relaxation time scale indicating the 

shortest period for processing all mass in the source layer is r sc - 0.33 h, but because TKEi is 
not sufficient to support the maximum mass flux, the actual shallow-convection period is -1.67 
h. However, in a weakly forced marine environment where the prevailing sensible heat flux is 
often < 10 Wm' 2 , it is common for TKE S < 1.0 J kg' 1 . In this situation, we may find D s -200 m 
and w B - 0.25 m s' 1 , so the relaxation time scale is t sc - 1 h, while the shallow-convection 
period becomes -10 h. Tests of the TKE-based closure showed it to operate well for most 
shallow clouds (D c < 2 km), while it damped most of the feedback oscillations that were 
characteristic of the BLD closure. 

2.3.5 Hybrid Closure 

Since the TKE-based closure was found to work well for fairly shallow convective clouds 
and the CAPE-removal closure was found to work better for deeper clouds (D c approaching 
D kf ), a simple hybrid closure is proposed to represent the intermediate range of cloud depths. 
When cloud tops are above the LFC, but have depths less than D KF (which describes a large 
percentage of convective cases), the clouds are assumed to be in transition from the TKE-based 
closure to the CAPE-removal closure. A simple linear averaging is used in this case, although 
refined transition functions could be hypothesized. First, the number of updrafts is calculated 
according to each of the two closures (N 2 and N 4 ). Then, the final number of updrafts is 
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estimated based on the fraction of the cloud depth that lies above the LFC relative to the distance 
between the LFC and the height D KF above cloud base, given by 

N = fN 2 +(l-f)N 4 (2.12) 

where f = /i 3 //i, is a ratio of two length scales (km), such that h 3 = D c -h 2 , h { = D KF -h 2 . 
Here, h 2 = z LFC ~ z B , and z LFC is the height of the LFC. Normally, N 2 <N 4 , because the 

CAPE-removal closure hypothesizes that stabilization of a deep cloud layer occurs rather rapidly 
as a result of a few vigorous updrafts. This agrees with the general observation that the number 
of growing clouds in an area decreases as their size (depth and width) increases due to 
compensating subsidence in the cloud environment. 

To summarize, the shallow-convection parameterization may use any of three mass-flux 
closure assumptions (type 2, 4 or 5), determined by D c , z T and z LFC , according to 

z T < z^c — » TKE - based closure only 

Zj > Zuc and D c < D KF —> Hybrid closure (2.13) 

D c > D kf — > CAPE - removal closure only 

which provides a smooth transition from one closure to another as the cloud depth grows. It also 
should be remembered that the uniform updraft geometry used in the present formulation 
represents a simplification, while the atmosphere typically exhibits an ensemble of cloud depths 
(see Sec. 2.2). Using LES, Siebesma and Cuijpers (1995) showed that a variety of cloud sizes 
contributes to a reduction in the updraft mass flux with height. While this effect is not 
considered in the present version of the parameterization, their LES results suggest that the 
updraft mass flux could perhaps be modified to account for a distribution of cloud sizes. 

2.4 Prognostic Scheme for Neutrally Buoyant Cloud Fraction and Cloud Water 

Most cloud parameterizations designed for the mesoscale consider only deep convection 
in detail. In these schemes, detrained cloudy air generally is fed back directly to the resolved 
scale, where it evaporates immediately until the grid-cell saturates (e.g., Kain and Fritsch 1993). 
This "all-or-nothing" approach for post-convective layer-cloud formation is an oversimplification 
and can have negative impacts on other aspects of model performance, such as radiative 
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processes and latent heating. In deep thunderstorms, large detrainment rates near the tropopause 
often cause a small to moderate-sized grid cell (DX < 25 km) to saturate quite rapidly, so that 
little damage is done by neglecting sub-grid layer clouds. However, in the case of smaller 
shallow clouds, it may take many hours (if ever) for detrained cloud water to saturate a layer. 

A more realistic representation is to detrain convective cloud mass from updrafts into a 
class of sub-grid clouds having nearly neutral buoyancy (e.g., Wyant et al. 1997). Once the 
detrained updraft air becomes part of these sub-grid NBCs, it can spread as layer clouds, initiate 
light precipitation, or slowly evaporate into the sub-saturated grid volume. While some existing 
schemes have attempted to treat this detrained cloud mass in large-scale models (e.g., Tiedtke 
1989, 1993), they often rely on a moisture-balance closure and so may not be versatile enough 
for both continental and marine environments. 

2.4.1 Basic Equations 

The mass-conserving equations for the rate of change of sub-grid cloud area (a) and 
cloud water/ice content (/ r ) for the NBCs are given by: 

^T = S a +D a ( 2 . 14 ) 

dt oz 

^ = s,+ D mix + D pre + D ics + D CTEl - V • V„/ f - ,| (2. 15) 

Here, S a and S, are sources of cloud area and condensed water content ejected from convective 
updrafts, respectively. The term D a (dissipation of cloud area) represents evaporation due to 
mixing at the sides of the cloud; D mix is the depletion of water liquid/ice content due to vertical 
mixing; D is water depletion due to precipitation (drizzle); D ics is a depletion rate contributed 
by an ice settling process and D CTEl is water depletion due to cloud-top entrainment instability. 
The horizontal and vertical advection terms are represented by - v • V H X an d ~ wdx/dz , where 
X is either l c or a . The grid-averaged condensed water content (/) is related to the sub-grid 
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NBC water content ( l c ) according to / = al c . When a = 1.0 (saturated conditions), then / = l c = 
q c , where q c is the model's explicit (resolved-scale) cloud-water mixing ratio. 

2.4.2 Formation of Neutrally Buoyant Clouds 

Convection produces a variety of clouds either directly, such as cumulus and 
cumulonimbus, or indirectly, such as stratocumulus and anvils. A realistic parameterization for 
clouds of convective origin, but which have nearly neutral buoyancy, is obtained by considering 
their source to be condensates produced in sub-grid cumulus updrafts (active convective clouds) 
and later detrained at the sub-grid scale into the non-convective environment. The following 
equations represent the source terms for the sub-grid NBCs described in (2.14) and (2.15): 

5 =-w— + (2.16) 

dz M i 

S,=(f-w/« (2.17) 

where ^- = -w-t-h — ~K (2.18) 

dt dz M l 

w is the convection-induced sub-grid-scale subsidence, R ud (kg [air] s '), is the updraft detrain- 
ment rate from the parcel-buoyancy scheme, and M L { kg) is total mass of air in a grid cell for a 
given model layer. The liquid/ice water content in the updraft core is given by l u . Note that the 
volume of detrained updraft air increases the area of the NBC (a), while its liquid water content 
( l c ) is solved as a residual term in (2.17) to satisfy the mass-conservation constraint l = al c . 

2.4.3 Evaporation of Clouds 

In (2.14) and (2.15), there are several processes through which the NBC can dissipate. 
Following Tiedtke (1993), the area decreases through cloud-edge evaporation, according to 

D„=-yK(q,-q,) (2.19) 
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where K = 10' 5 s' 1 is a horizontal diffusion coefficient; and q v and q s are the resolved-scale 
specific humidity and saturation specific humidity, respectively. This yields a dissipation time 
scale on order of a day for l c = 1.0 g kg' 1 air. The dissipation rate becomes significantly greater, 

however, when the cloud water content becomes small (l c < 0.1 g kg' 1 air). Thus, in effect, this 
term accounts for the final dissipation of clouds from which most of the condensate has already 
been depleted by other processes. 

2.4.4 Depletion of Water Liquid/Ice Content by In-Cloud Mixing Processes 


Cloud water in the sub-model can be mixed downward (but not upward) inside the NBC 
by vertical diffusion induced by turbulence and radiation flux divergence, given by 


n ^ ^ UK i K 1 


( 2 . 20 ) 


where K is the local diffusion coefficient derived from the TKE. The additional radiation- 


induced diffusion coefficient is given by 


K =i\^- ^JL+— — 1 


0 L dt 


15 dt sw 50 


( 2 . 21 ) 


where K r is a maximum at the cloud top and decreases linearly downward over a maximum 
cloud depth of 1000 m. In (2.21), l B is the Blackadar length scale provided by the turbulence 
scheme (Shafran et al. 2000), 0 is the potential temperature of the environmental air and A; r is 

dO . d0 . . 

the model layer thickness at cloud the top. The terms — and are the longwave 

Ol r w or w 


cooling rate at cloud top and the daytime solar heating just below cloud top, respectively, 
provided by the radiation scheme. This K r term would be unnecessary if very high vertical 
resolutions were possible, but in typical mesoscale models with a Zt ^ 50 m, it accounts for 
radiation-induced turbulence occurring at scales that too small to be represented accurately by 
K v , which is derived from grid-resolved quantities in the turbulence scheme. 
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Cloud water also can evaporate into clear air at the exposed cloud base through the same 
vertical diffusion terms. The finite-differenced form of these diffusion terms is one-sided so that 
the water flux is downward only, while mixing into layers above the cloud is prohibited. These 
vertical mixing terms become the most important water depletion process when the liquid/ice 

content is in the range of 0. 1 < l c < 0.5 g kg" 1 , where neither sub-grid precipitation (Sec. 2.4.5) 
or horizontal diffusion processes (Sec. 2.4.3) are effective. 

2.4.5 Precipitation Processes 

Precipitation can form in the NBCs through auto-conversion and accretion, exactly as for 
resolved-scale layer clouds. Here, we use the simple water/ice cloud scheme of Dudhia (1989), 
without a mixed phase, for which auto-conversion begins at l c >0.5g kg' 1 (a similar threshold 

exists for initiating auto-conversion of ice, based on activation of ice nuclei below 273 K). It is 
feasible to introduce a mixed-phase explicit precipitation scheme for use with the shallow-cloud 
scheme, but that is not done here. Also, cloud-ice particles are assumed to fall slowly through an 
ice-settling process described by Grell et al. (1994). The ice settling prevents long-term retention 
of low-density cirrus clouds that are unable to produce significant snowfall rates. 

2.4.6 Cloud-Top Entrainment Instability (CTEI) 

Deardorff (1980) proposed that shallow clouds can dissipate through a mechanism called 
cloud-top entrainment instability (CTEI). If a parcel of dry air is entrained into the cloud top, it 
induces mixing and evaporation. As a result, the density of the parcel may become greater than 
(G v less than) that of surrounding cloudy air, causing unstable acceleration of the parcel 
downward through the cloud. The depletion rate of the cloud liquid due to the CTEI, following 
Deardorff, is parameterized as 

£crc/=-10~ 4 — ° — l c (2.22) 

r — V 
'max 'min 
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Deardorff (1980) shows that r = ^ h / L ^ q + / y where the moist static energy ,S glV6n by 
/j = Q p j + g Z + Lq v , A h represents the jump of h across the cloud top, A (q v +l c ) is the jump 
of total water across the cloud-top, L is the latent heat of vaporization, and C P is the specific 
heat at constant pressure. In (2.22) and r max are CTEI initiation criteria given by Deardorff 
(1980) and MacVean and Mason (1990). If r < r mm , then r is set equal to r IIun and the depletion 
rate goes to zero; if r>r max then r is set equal to r max , which gives the maximum depletion 
rate, -10~ 3 4 / c kg kg'V . Because CTEI-induced downdraft cooling is expected to be confined 
to the upper portion of the cloud layer (Randall 1980), the depletion due to CTEI is applied only 
in the uppermost 100 m of the NBC. 

Deardorff (1980) originally proposed CTEI as a possible dominant mechanism tor the 
breakup and evaporation of a stratus deck. Although subsequent observational and modeling 
research has suggested that this is generally not the case, the process is included in the NBC sub- 
model as a contributing factor for cloud water depletion. 


3. EXPERIMENTATION IN A SINGLE COLUMN MODEL 
3.1 The 1-D MM5 

The initial testbed for the shallow-convection parameterization has been a 1-D version of 
the Penn State/NCAR mesoscale model (MM5), which is based on the non-hydrostatic numerical 
framework of the 3-D parent model (Grell et al. 1994). In addition to the shallow-convection 
scheme, the 1-D model contains a full set of MM5 physical parameterizations, including a 1.5- 
order turbulence sub-model (Shafran et al. 2000), an explicit moisture scheme (Dudhia 1989) 
that predicts resolved-scale liquid/ice cloud and precipitation, and a two-stream broadband 
column radiation sub-model (Dudhia 1989). The land-surface scheme can be chosen to be either 
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the Blackadar Slab model or PLACE. As discussed in Sec. 2, it also contains the KF (1990) 
deep convection parameterization. Moreover, the user can specify 3-D large-scale dynamical 
tendencies, such as advection and subsidence. Of these, the interaction between the shallow 
convection and radiation requires special note. 

To function correctly with the shallow convection sub-model, the radiation scheme must 
be applied separately to the sub-grid clear and cloudy areas of the column. Like many other 
schemes, however, Dudhia's radiation assumes that a grid cell is either totally clear or totally 
cloudy. Furthermore, the shallow-cloud area predicted by (2.14) is a function of height. To 
solve this problem efficiently, the column's sub-grid NBC is partitioned into three idealized parts 
(see Fig. 3.1.1): (1) clear throughout the model depth, (2) sub-grid upright cloud through the full 
depth of the convective layer (often this is a fairly small fraction of the grid area), and (3) broad 
stratiform sub-grid cloud (often near the top of the convective updraft, referred to here as the 
“anvil” for convenience). The radiation-induced thermal tendencies are calculated and applied 
separately to each fraction of the cell area. The convective updraft is detrained completely into 
the NBC after each time step, and so does not contribute to the radiation tendencies. Deng 
(1999) gives details of the methodology used for the sub-grid partitioning represented in Figure 
3.1.1. 

If the partitioned sub-grid cloud areas 2 and 3 in the figure (a 2 and a 3 ) were used 

directly for the radiation calculations, however, the true cloud effects would be underestimated. 
The "effective" cloud fraction (i.e., the area as viewed from above) must be somewhat greater 
than the idealized distribution of Fig. 3.1.1 because (2.14) implicitly assumes that clouds in 
adjacent layers always have the maximum possible overlap in the vertical. While this maximum- 
overlap assumption gives a reasonable first guess for clouds of convective origin, it is not likely 
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Figure 3.1.1. 


Schematic diagram showing idealized partition of NBC in a grid column into (1) 
clear, (2) upright cloud, and (3) anvil cloud. 


35 




to be accurate because in nature many dissipating NBCs are no longer linked vertically by an 
active updraft. An additional fractional area, a s , is postulated to be a function of the 

environmental relative humidity, RH , so that the total or "effective" cloud fraction, a e , can be 
written as 

a e — (1 — u)a s + a (3-1) 

where a is the fraction predicted using (2.14). Following Xu and Randall (1996), the additional 
fraction, a s , is parameterized according to 

= RH a ' [1 - exp( ~ a ^ lc+q ^ )) if RH < 1 (3.2) 

s F [(l -RH)q s p 

where a s = 1 if RH > 1 and or, , a 2 and <2 3 are empirical constants defined as 0.25, 0.49 and 
100, respectively. Note that the cloud liquid/ice content used in (3.2) appears as the sum of the 
cloud liquid/ice at the sub-grid scale, l c , and the resolved scale, q c . However, when a layer 
becomes saturated, l c converts into q c , while q c — 0 when the layer is sub-saturated, so only 

one of these two quantities can actually be non-zero at a given time. 

As a brief example, Fig. 3.1.2 shows the resultant area of effective cloud, a e (thin solid 
curve), in a marine stratocumulus case where the maximum predicted cloud area from (2.14) is 
fl max =0-44 (shaded region) and RH m ax =0.93 (heavy solid curve) at 1451 m. Notice that, in 
the cloud layers with lower relative humidity (-0.85-0.90), the effective cloud fraction is only 
about 5-10% greater than the calculated NBC fraction. However, near the cloud top, where the 
relative humidity is greatest, a e exceeds the NBC fraction by nearly 40%. 

The surface radiation flux, R s , (for either LW or SW flux) is given by 

R$ — ^ e 2^2 4" ^e3^3 (1 — @e2 ~ (3-3) 
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Figure 3.1.2. Model-predicted cloud fraction (%) and resolved-scale relative humidity, RH (%), 
versus pressure (mb) and height (m) at 1200 UTC, 11 June 1992 (Hour 5) for 
ASTEX at 28.00°N, 24.22°W. Shading indicates NBC fraction (%), the dashed 
line is the area of shallow-cloud updraft (%), the thin solid curve is the effective 
cloud area for radiation calculations (%), and the heavy curve is RH (%). LCL 
denotes the lifting condensation level, PBL is the boundary layer top, CLDTOP is 
top of updraft, and NBCT and NBCB are the top and base of NBCs, respectively. 
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where a e2 and a e3 are the vertically averaged effective cloud fractions and R 2 and R 3 are the 
surface fluxes contributed from the sub-grid upright and anvil portions of the NBC, respectively. 
The term R clr is the surface flux contributed by the clear portion of the grid element. Finally, 
Fig. 3.1.3 presents the evolution of the surface radiation fluxes for each of these three sub-grid 
components and the net grid-averaged surface radiation determined by (3.3) for three days during 
the same marine stratocumulus case (11-14 June 1992) from which Figure 3.1.2 is taken. 

Naturally, calling the Dudhia (1989) radiation sub-model three times for each grid cell to 
account for sub-grid cloud effects involves extra calculations, compared to the usual approach 
when there are no sub-grid clouds. The added computational load due to the multiple 
applications of the radiation module can be offset in part by not calling the scheme each time 
step. It should be sufficient to call it about once every 10-15 minutes, unless the cloud field is 
changing very rapidly. Furthermore, multiple calls to the radiation module could be avoided, in 
principle, by developing an equivalent grid-averaged representation for the radiative properties 
of the sub-grid cloudy and clear portions of the column. However, this latter concept is more 
complex than merely spreading the sub-grid cloud water horizontally across the grid cell. Such 
an approach is beyond the scope of the present study. 
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TIME (UTC) 


Figure 3.1.3. Model surface solar fluxes for 0600 UTC, 1 1 June to 0600 UTC, 14 June 1992 for 
ASTEX at 28.00°N, 24.22°W. Thin solid line is the clear-sky solar flux ( R clr ), 

dashed line is the solar flux for upright cloud ( R 2 ), dotted line is solar flux for the 
anvil cloud (/?,), and heavy solid line (R s ) is the area-weighted average flux. 




3.2 Offline (Stand-Alone) Experiments with PLACE 


3.2.1 Introductory Remarks 

We describe in this section the offline (stand-alone) PLACE-related work for this study. 
The offline version of the PLACE model was forced by meteorological observations (rather than 
an atmospheric model) to create the initial conditions for the 3-D coupled MM5-PLACE 
simulations discussed in Section 4. The objectives of this work were (1) to produce realistic and 
consistent initial fields for the soil variables required by the coupled 3-D MM5-PLACE model, 
and (2) to examine the overall performance of PLACE over a large domain and over a long 
period when forced with analyzed data and using currently available characterizations of soil and 
vegetation fields. This offline PLACE modeling system and its results are described below. 

3.2.2 Model Description 

The stand-alone version of PLACE was constructed starting from the real-time Soil 
Hydrology Model (SHM) modeling framework already available at Penn State (Smith et al. 
1994). The SHM system produced an archive of analyzed meteorological fields that had been 
used to feed the SHM model in real time for several years, including the summer of 1997, which 
is the period that we have focused on for testing and application of the MM5-PLACE model. 

The available fields provided by the SHM archive are 12-h precipitation, daily cloudiness, daily 
values of wind speed and indirect measures of temperature and humidity. In order to produce a 
more realistic representation of the diurnal cycle in the forcing, the SHM-preprocessor programs 
were re-run to construct meteorological forcings every three hours for wind speed, air 
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temperature, water mixing ratio and surface pressure. A simple radiative model, also from the 
SHM, is used to calculate radiative forcings for PLACE based on the cloud fields. 

Figure 3.2.1 shows a schematic of the PLACE stand-alone model. Whenever possible the 
surface characterization fields and parameters used in the stand-alone runs were the same as 
those available in the standard MM5V3 modeling system, so that the results could be used 
directly as initial conditions in the coupled MM5-PLACE model runs. The next subsection 
compares the analyzed forcing fields with available measurements in the main area of interest. 

3.2.3 Verification of the Forcings 

The most important forcing for a land-surface model is precipitation. In Fig. 3.2.2a, we 
compare for the month of July 1997 the time series of daily precipitation amounts averaged over 
110 Oklahoma Mesonet Stations with the analyzed precipitation amounts coming from the SHM- 
real-time system, for the grid points closest to the station locations. (The Mesonet data were 
acquired through another related NASA grant to study the effects of surface heterogeneity during 
the Southern Great Plains 1997 [SGP97] study period). Conventional analyses were not 
available for days 14, 15, and 22, so fields from the previous days were repeated to fill the gaps. 
The main precipitation events over the region are well captured in the analyzed data. An 
apparent exception of this is the precipitation amounts for 9 July 1997. The analysis shows 
significant precipitation that was not measured at the stations. Study of the horizontal 
distribution of the analyzed precipitation and of radar data shows that indeed there existed a 
significant precipitation event close to the limits of the Mesonet region, so that a small location 
error in the analysis appears to have produced a large bias in the total precipitation of the 
analysis for the Mesonet points. Figure 3.2.2b shows a scatter plot of daily averaged 
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Figure 3.2.1 Schematic of PLACE Stand-Alone Model 
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Figure 3.2.2 Comparison of forcings of stand-alone PLACE model and Mesonet station 
measurements for July 1997. a) Time series of measured and analyzed mean daily precipitation 
amounts, b) Scatter plot of mean daily precipitation for all stations, c) Time series of measured 
and analyzed mean daily solar radiation, d) Scatter plot of mean daily solar radiation for all 
stations, e) Time series of measured and analyzed mean daily longwave radiation at station E13, 
f) Time series of measured and analyzed mean atmospheric surface pressure. 
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precipitation measured and analyzed for the same 110 Mesonet stations with precipitation data. 
The daily averages are computed for the full month of July 1997. Although the figure shows 
considerable scatter, there is no significant bias between measured and analyzed precipitation. 

Figures 3.2.2c and 3.2.2d show a similar analysis, but for analyzed and measured solar 
radiation (67 stations with full solar radiation measurements are used in this analysis). The 
analyzed solar radiation data appear to underestimate the shortwave radiation reaching the 
ground for clear conditions and overestimate it in cloudy conditions. Figure 3.2.2e is a check for 
the longwave radiation reaching the surface. Longwave radiation is not measured over the 
Mesonet stations, so that the comparison shown in the figure is restricted to station E13 of the 
ARM-CART site (siros datastream). The model longwave radiation reproduces reasonably well 
the measured flux. Finally, Fig. 3.2.2f compares time series of atmospheric surface pressure 
averaged for the Mesonet stations and the corresponding analyzed fields. The mean surface 
pressure bias is around 2-3 hPa. These analysis errors may be related in part to the relatively 
coarse 36-km terrain fields in the offline PLACE domain (Fig. 3.2.3), compared to the actual 
terrain at the Mesonet sites. 

Figure 3.2.4 shows temporal and spatial verifications of meteorological forcings of the 
PLACE stand-alone system. Analyzed fields of air temperature (panels a and b), water vapor 
mixing ratio (panels c and d), and wind speed (panels e and f) agree reasonably well with the 
independent high-resolution measurements available in the Mesonet. Although the forcings that 
feed the offline PLACE may be further improved by using the Oklahoma Mesonet data, the 
offline PLACE must be run over the entire 36-km continental-scale MM5 domain (see Fig. 

3.2.3) where special data do not exist, and these Mesonet data will also serve as an independent 


44 






Mean Temperature Time Series 





30 


29 


28 


■ 27 


26 


I 25 
w 

24 

24 


Observed vs Analyzed Temperature 


b 









• , <• 

* : 

■ 








: 


















25 26 27 28 29 

Mesonet Temperature (C) 


30 


Observed vs Analyzed Mixing Ratio 



Mesonet Mixing Ratio (g/g) 


Observed vs Analyzed Wind Speed 



0 1 2 3 4 5 6 

Mesonet Wind Speed (m/s) 


Figure 3.2.4 Comparisons of forcings to stand-alone PLACE model and stations measurements 
for July 1997. a) Time series of measured and analyzed mean temperature, b) Scatter plot of 
mean daily temperature for all stations, c) Time series of measured and analyzed mean water 
vapor mixing ratio, d) Scatter plot of mean daily water vapor mixing ratio for all stations, e) 

Time series of measured and analyzed mean wind speed, f) Scatter plot of mean daily wind 
speed for all stations. 
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verification of the 3-D MM5-PLACE surface fields in section 4.1. Thus, this independent 
verification of the conventional analyses used to force offline PLACE against the Mesonet data 
shows that these data are reasonably accurate for the purpose of creating the initial conditions for 
the MM5-PLACE 3-D experiments. 

3.2.4 Model Settings and Parameter Specifications 

The model domain (Fig. 3.2.3) is based on the SHM real-time system (Smith et al. 1994). 
It covers most of the continental United States with 94 points in the N-S direction and 1 14 points 
in the W-E direction. All the work presented in this report uses this 36-km domain. The time 
step used in the stand-alone runs of PLACE is 10 minutes. Although PLACE has the capability 
of handling several mosaic tiles in each grid cell, only one tile has been used in these runs. 
Whenever possible the surface characterization and the soil and vegetation parameters are taken 
from the MM5V3 system. The soil is characterized based on the 16-category STATSGO soil 
map. The land-use/vegetation characterization is based on the 24-type USGS classification 
scheme that is available as a standard MM5 option. Some of the parameters of PLACE that are 
not specified in the MM5V3 land-surface model have been defined with the help of the PLACE 
developer, Dr. Peter Wetzel. Among these parameters are the vegetation leaf area index (LAI) 
that has been set at a value of 7. The map of vegetation fraction used in these runs comes from 
the climatological vegetation fraction map available in MM5V3 for the month of July. A more 
detailed description of the parameters required by PLACE, and their definition in the results 
presented here can be found in Appendix A and Munoz (2002). The results presented in 
subsequent sections correspond to an offline model run that was initialized on 1 June 1997. 

The PLACE model uses 7 soil temperature layers defined at 0-2 cm below the surface, 2- 
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5 cm, 5-10 cm, 10-15 cm, 15-50 cm, 50-90 cm and 90-130 cm. Its 5 soil moisture layers are 
located at 0-2 cm, 2-5 cm, 5-15 cm, 15-50 cm and 50-130 cm. The initial soil moisture field was 
derived from the moisture field at 1 June 1997 of the SHM-real-time-model. The soil moisture 
profile at the PLACE levels was obtained by linear interpolation of the SHM soil moisture 
profile. The initial soil temperature profile was calculated as a linear interpolation between a 
surface soil temperature and a deep soil temperature field. The upper field is calculated as the 
daily average of the air temperature, and the deep soil temperature is the climatological deep soil 
temperature field available in the MM5V3 system. 

3.2.5 Offline PLACE Results 

Results presented in this section are obtained using the model settings and parameters 
described above, with the only exception that the precipitation forcing used was obtained from 
the measurements at each individual Mesonet station. This procedure allowed us to perform a 
more extensive comparison between model results and point measurements than if we used the 
analyzed precipitation fields. Although the latter fields had not shown any significant bias with 
respect to the measurements (except July 9, as noted earlier), the variance between 
measurements and analysis made the evaluation of the PLACE model at the station level less 
revealing. However, the gridded fields of soil variables used as initial conditions for the 3-D 
MM5-PLACE runs are produced by using the analyzed precipitation fields as forcing for the 
stand-alone PLACE model. 

3.2.5. 1 Soil Moisture 

Figure 3.2.5 shows daily time series of observed and modeled soil moisture averaged 
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THETA05 Daily Values 1-30 July 97 


THETA25 Daily Values 1-30 July 97 






Figure 3.2.5 Time series of measured and modeled volumetric soil moisture content averaged 
over 35 Mesonet stations with full soil moisture data for July 1997. a) 5 cm depth, b) 25 cm 
depth, c) 60 cm depth, d) 75 cm depth. 
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over 35 Mesonet stations, for the period 1 to 30 July 1997, and for the four levels with 
measurements (model results for each station are shown in Appendix B). Soil moisture at 5 cm 
depth (panel a) in the measurements and the model results responds strongly to the precipitation 
events during the month. The rain episodes produce increments of similar magnitude in the 
model as in the observations. The drying trend in between rain events, however, is different in 
the model as compared to the data. The data show larger drying rates at this level. The second 
level at 25-cm depth (panel b) shows a remarkable similarity between the observed and measured 
mean drying trends during the month. Both time series averaged over the 35 Mesonet stations 
show a systematic difference that is due to an initial error and to the effects of the rain event 
around July 10 that is significant in the model results, but is only slightly noted in the 
observations. Mean results for the 60-cm deep soil moisture (panel c) are also similar between 
model and observations, again with a slightly larger drying in the observations than in the model 
results. The same is evident at the depth of 75 cm (panel d). 

Figure 3.2.6 shows number histograms of soil moisture at different levels, as derived 
from model results and observations. The difference in the shapes of the histograms is 
significant. Model results show a tendency to peak at middle values of soil moisture, while 
observed histograms are highly skewed with larger frequencies at the low soil moisture range. A 
secondary maximum in the observed histograms occurs at the high soil moisture range, 
especially in the deeper soil layers. Therefore, the PLACE model is having some difficulty 
simulating the very low and very high soil moisture values, especially in the deeper soil layers. 
Deep soil values may be improved by extending the length of the PLACE simulation to several 
months. 
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Figure 3.2.6 Number histograms of measured and modeled soil moisture contents for July 
1997. a) 5 cm depth, b) 25 cm depth, c) 60 cm depth, d) 75 cm depth. 
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Figure 3.2.7 Scatter plot of monthly mean soil moisture for 35 Mesonet stations with full soil 
moisture data for July 1997. a) 5 cm depth, b) 25 cm depth, c) 60 cm depth, d) 75 cm depth. 






Figure 3.2.7 shows the dispersion between monthly mean soil-moisture values from the 
measurements and the model results, for the four levels with data. It is apparent that the model is 
unable to capture the regions with larger soil moisture values in the data. This bias could be due 
to excessively rapid drying in the model or simply due to soil moisture parameters that are 
different between model and observation sites (for example, a saturation moisture that is smaller 
in the model than in the measurements). We investigate this problem further by displaying the 
dispersion between the observed and modeled maximum and minimum soil moisture at each 
station, as shown in Fig. 3.2.8. Panel a shows the scatter between monthly minimum soil 
moisture at 5-cm depth (minimum values are paired independent of the day that they have been 
reached). The minimum values of soil moisture are in many cases significantly greater than the 
corresponding values in the observations. The scatter of maximum values at 5-cm depth (panel 
b) shows a smaller bias between model and observations in the middle soil moisture range. The 
largest soil moisture amounts, however, are consistently underestimated by the model. Similar 
conclusions are drawn for the level 25-cm deep, as shown in panels c and d. 

Figure 3.2.9 shows time series of the correlation coefficients between observed and 
modeled soil moisture for the four depths with data. A significant increase in the model- 
observation correlations is evident during the rainy episodes, especially in the two uppermost 
levels. The two lower levels show a steady increase in the correlation through the monthly 
period, suggesting that a longer model run could improve the correlations between model and 
measurements at the deeper levels. The two lines in these panels correspond to the sensitivity to 
the precipitation forcing used in the model. Run 14s8 (continuous line) used the observed 
precipitation rates starting on 1 July 1997 (during June 1997 it used the analyzed precipitation 
fields), while run 14s8sl used the measured precipitation rates from the start on 1 June 1997. 
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Figure 3.2.8 Scatter plot of monthly minimum and maximum soil moisture for 35 Mesonet 
stations with full soil moisture data for July 1997. a) minimum soil moisture at 5 cm depth, b) 
maximum soil moisture at 5 cm depth, c) minimum soil moisture at 25 cm depth, d) maximum 
soil moisture at 25 cm depth. 
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Figure 3.2.9 Time series of correlation coefficients between observed and measured soil 
moisture contents for the 35 Mesonet stations with full soil moisture data for July 1997. Model 
run Runl4s8 is driven by point-measured precipitation starting on 1 July, and run Runl4s8sl 
uses point-measured precipitation starting from 1 June, a) 5 cm depth, b) 25 cm depth, c) 60 cm 
depth, d) 75 cm depth. 
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The effect of using better precipitation forcing is evident in the lower layers. The 
uppermost layer shows larger correlation coefficients only during the first days of July. The 
monthly correlation between model and observations for each station is represented in Fig. 

3.2.10 for all four levels. The uppermost levels show only positive correlations. The lower 
levels show more variable correlations, with some stations showing negative values and others 
showing values closer to 1. 

Another way of studying the spatial variability in the observations and in the model 
results is by grouping the stations in four regions, each covering (overlapping) regions of 
approximately 2x2 degrees. These regions are defined as A, B, C, D, from the westernmost to the 
easternmost part of Oklahoma. Figure 3.2.1 1 shows the mean time series of soil moisture in the 
four regions according to the available measurements. A significant west-east gradient in soil 
moisture is apparent from this figure. The western portion of the Mesonet is significantly drier 
than the eastern part. This gradient is more noticeable in the second layer (panel b). The upper 
layer is more affected by the rain distribution, and the lower layers appear to be more 
homogeneous in soil moisture, with the exception of region A that is drier in all four levels. The 
corresponding regional averages from the model results are presented in Fig. 3.2.12. The west- 
east gradient is not well reproduced in the model results. Only at the end of the month has a 
gradient developed after the large precipitation event around 19 July 1997. 

The mean vertical soil-moisture profile is shown in Fig. 3.2.13, where the mean profiles 
are plotted for observations and model results. A mean overestimation in the soil moisture 
profile is evident from the figure, although the vertical gradients of soil moisture are similar 
below 20-cm depth. The model results show a systematic minimum in the second level. 
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Figure 3.2.10 Monthly correlation coefficients between model and measured soil moisture 
contents for 35 Mesonet stations with full soil moisture data for July 1997 (station index does not 
correspond to Mesonet station number), a) 5 cm depth, b) 25 cm depth, c) 60 cm depth, d) 75 cm 
depth. 
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Figure 3.2.1 1 Time series of measured volumetric soil moisture content at 35 Mesonet stations 
with full soil moisture data for July 1997 and for the four regions described in the text (Region A 
is the westernmost and Region D is easternmost), a) 5 cm depth, b) 25 cm depth, c) 60 cm depth, 
d) 75 cm depth. 



Days July 97 Days July 97 


Figure 3.2.12 Same as in Fig. 3.2.10, but for model results. 
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This may be due to a possible specification error in the root fraction at this level (0.5), which 
may make it too large. Further experimentation with longer PLACE integrations is necessary to 
obtain definitive conclusions. 

Finally, it has been possible to compare the offline PLACE model results to high 
resolution surface soil moisture derived from the Electronically Scanned Thinned Array 
Radiometer (ESTAR) instrument flown on a P3 aircraft during SGP97 (Jackson et al., 1999). An 
example of this comparison is shown in Fig. 3.2.14 taken from Reen et al. (2001). The figure 
shows a cross section of soil moisture content along one of the ESTAR southwest-northeast 
flight paths on 12 July 1997. The cross section is around 400 km long and the ARM-CART 
Central Facility at Lamont, Oklahoma, is located near grid cell 95 in the figure. The 800-m 
resolution ESTAR data and the offline PLACE data were both averaged to a 4-km grid for this 
comparison. On this day there was a general north-south gradient in soil moisture with warmer, 
drier conditions towards the south largely due to the convective precipitation that occurred the 
day before in the northern region. In spite of the coarse 36-km resolution of the PLACE stand- 
alone results, the spatial variability at this scale is remarkably well captured by the model results 
on this day. 


3.2.5.2 Soil Temperature 

Most Mesonet stations have measurements of soil temperature at 10-cm depth and about 
half of these also have measurements at 5-cm and 30-cm depths. Panels a, b, and c of Fig. 3.2.15 
show monthly time series of the soil temperature at different levels averaged over all stations 
with data (and over all 111 Mesonet stations in the case of model results). The two types of data 
lines correspond to measurements below bare soil (TB05 and TB10) and below a sod 
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Figure 3.2.14 Offline PLACE versus ESTAR soil moisture contents for 12 July 1997. Each 
model cell y-value represents an east-west 4-km grid average across the width of the ESTAR 
data area (from Reen et al. 2001). 


61 




Soil temperature at 5 cm Soil temperature at 10 cm Soil temperature at 30 cm 




Diurnal cycles of soil temperature at DSucnal cycles of soil temperature at DQUcnal cycles of soil temperature at 30 cm 

40 r 




35 


30 


25 


20 


f 

-4- TS30 
R14s8 



44 t F4 i i i 

-ri • ■ ^ ■ i m i i i | | | | | | | | | | r 


0 6 12 18 24 

Hours 


Figure 3.2. 15 Comparison of observed and modeled soil temperature, a) Soil temperature at 5 
cm depth averaged over all Mesonet stations with data and all 111 stations in the case of model 
results. TB05 is measured below bare soil and TS05 is measured below a sod cover, b) as a) but 
for 10 cm depth, c) as a) but for 30 cm depth, d) Mean diurnal cycles for July 1997 at 5 cm 
depth, e) as d) but for 10 cm depth, f) as d) but for 30 cm depth. 
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cover (TS05, TS10, and TS30). The measured mean temperatures below bare soil are 
consistently higher than below sod cover, also have a larger diumal amplitude (not shown). The 
model soil temperatures in the upper two layers fall between the two measured averaged series, 
with the exception of the period around 6 July, when the model temperatures dropped more than 
the observations. The slight warming trend of the observed deep soil temperature is slightly 
smaller than the monthly trend of the model results. Panels d, e, and f of the same figure show 
diumal cycles of measurements and model results at the three levels with data. The amplitude of 
the mean diumal cycles in the model agrees well with the measured cycles at the different soil 
depths. 

3.2.5.3 Surface Fluxes 

Comparison of modeled and observed surface fluxes is always complicated by the local 
character of the measurements. Mesonet stations do not have flux measurements. Stations in the 
SGP ARM-CART site include EBBR and eddy correlation measurements of fluxes. Many 
stations, however, present measurement problems during extended periods. We restrict our 
comparison here to ARM-CART stations E9 and E15 that do not report data quality problems 
during July 1997. Figures 3.2. 16 and 3.2. 17 summarize the analysis of the sensible and latent 
heat fluxes measured at these stations. Diumal cycles of fluxes are shown in panels a and c of 
Fig. 3.2.16. In both stations the surface fluxes appear to be dominated by the latent heat fluxes, 
especially at station E9. Panels b and d illustrate the diumal cycle of the variability in these 
measurements. Noise in the data is relatively large during the night and in the diumal transition 
periods, where the Bowen ratio method of determining fluxes usually has problems. The fluxes 
between hours 15 and 21 are averaged for each day and presented in the monthly time series in 
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Surface Fluxes E9 



Standard Deviation of Fluxes E9 



Surface Fluxes El 5 



Standard Deviation of Fluxes El 5 



Figure 3.2.16 Measured surface fluxes at two ARM-CART EBBR stations, a) Mean diurnal 
cycles of sensible and latent heat fluxes for station E9 during July 1997. b) Mean diurnal cycles 
of standard deviation of surface fluxes for station E9 during July 1997. c) As a) but for station 
E15. d) As b) but for station E15. 
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Figure 3.2.17 Monthly time series of diurnal surface fluxes at two ARM-CART EBBR stations. 
Each daily value is the average of measurements between 15 UTC and 21 UTC. a) Station E9, b) 
Station El 5. 
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Figure 3.2.17. The dominance of latent heat fluxes at station E9 is again evident. The low 
values of sensible heat flux during the afternoon at E9 for the entire month is intriguing. There 
appears to be some source of water nearby. Station E15 also has Bowen ratios less than one for 
more than two-thirds of the month, but it has Bowen ratio near one (i.e., comparable latent and 
sensible heat fluxes) during the period of 5 to 15 of July 1997. 

The corresponding diurnal cycles and monthly time series of fluxes for offline PLACE 
model results are shown in Fig. 3.2.18. Mean model results for station E9 are also dominated by 
the latent flux, although not as strongly as in the observations. Station E15 has more comparable 
fluxes. The 3-hourly noise evident in these diurnal cycles is probably an artifact of the 
meteorological forcing supplied to the model that changes every three hours. Panels c and d of 
the same figure roughly follow the corresponding observed time series. The most evident 
problem in the model results is when the sensible heat flux remains very low and even negative 
during daytime. This may be due to the use of one cloud field for the entire day in the radiation 
model. The corresponding observed fluxes become small, but usually remain positive. Overall, 
the offline PLACE results are representing fairly well the land-surface conditions and its daily 
variability. 
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Surface Fluxed E9 Run14s8arm 



Surface Fluxed El 5 Run14s8arm 



Hours 15-21 Averaged Surface Fluxed E9 Run14s8arm 




Figure 3.2.18 Analysis of modeled surface fluxes for two ARM-CART EBBR stations, a) Mean 
diurnal cycle of surface fluxes at station E9 for July 1997. b) As a) but for station E15. c) As 
Figure 3.2.17a but for model results, d) As Figure 3.2.17b but for model results. 
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3.3 The Improved Shallow-Convection Scheme 


3.3.1 Removal of Stochastic Behavior in the Mass Flux Calculations 
3.3.1. 1 Introductory Remarks 

In order to use the shallow convection parameterization to investigate the interaction 
among shallow clouds and other physical processes like radiation, turbulence, deep convection 
and surface processes, it is necessary first to make sure that the parameterization is working in a 
physically plausible way. This means that, among other things, it appropriately conserves 
quantities like water or energy and that it produces realistic tendencies. It also means that it has 
a smooth behavior, in the sense that for small changes in the variables input to the scheme, it 
produces smoothly varying tendencies. 

Preliminary application of the shallow scheme has shown that the latter conditions are not 
fully satisfied and that several improvements to the shallow scheme are necessary if one wants to 
use it with the purposes of regional climate and atmospheric water-budget studies. This section 
describes the most important problems encountered in the original scheme of Deng (1999) and 
points to suggested solutions and improvements. 

For the purpose of this report we divide the shallow-cloud scheme into three main parts. 
Part 1 will be referred to as the updraft algorithm. It defines the initial parcel that is used to 
determine the convective amount (intensity) and calculates the updraft mass fluxes for that parcel 
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as a function of the given environmental sounding. Part 2 of the shallow parameterization will 
be referred to as the closure algorithm. It computes the environmental tendencies derived from 
the updraft mass flux profile and evaluates the change in available buoyant energy due to those 
tendencies. It also determines how much convection will occur, given an appropriate closure 
condition. The final part of the shallow-cloud scheme, Part 3, calculates the dissipation of 
shallow cloud variables in the neutrally buoyant cloud mass that is detrained from the shallow- 
cloud updrafts. The following portions of this section describe the problems and improvements 
for each of the three parts of the scheme and illustrate results using the original formulation 
versus the improved scheme. 


3.3.1.2 Improvements to the Kain-Fritsch Updraft Algorithm 

Preliminary runs with the shallow convection scheme have shown a large degree of 
discontinuous behavior in its solutions. This feature hinders the attainment of quasi-equilibrium 
conditions, as well as making the physical interpretation of model results quite difficult. One of 
the underlying reasons for the discontinuous behavior of the shallow convection scheme is that 
its cloud-updraft model is based on the Kain-Fritsch (1990) deep-convection scheme. The latter 
was designed to calculate a mean updraft condition to be used for a convective period lasting on 
the order of a half hour or so. Therefore, the discontinuous features of its cloud model had little 
impact over a period of several time steps. When a new Kain-Fritsch updraft was determined a 
half hour later, it would be calculated on the basis of what could be (by that time) a very different 
environment. However, the original shallow-convection scheme calculates a new updraft at 
every time step using this same Kain-Fritsch cloud model, which causes the discontinuous 
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tendencies that become painfully evident. 


3.3.1.2.1 Case 0: A Simple Case with Discontinuities 

As the first step in the investigation, we isolated the algorithm that computes the source 
conditions for the shallow clouds. In essence it corresponds to the source algorithm in the Kain- 
Fritsch (1990) cloud model, with some modifications regarding the source layer and the initial 
level of the updraft parcels defined for the shallow convection scheme. To test the scheme we 
used the thermodynamic sounding shown in Figure 3.3.1. It consists of a 500-m mixed layer 
capped by a deep isothermal layer. The water vapor mixing ratio is constant in almost all the 
lower atmosphere above 990 mb. For this case the lowest two layers of the sounding, 
comprising a depth of about 10 mb (-100 m), define the thermodynamic properties of the updraft 
parcel. The mixing ratio of these two layers is moister than in the upper levels, but is kept 
constant in all calculations shown in this section. The temperature of the two lowest layers will 
be varied to study the continuity (smoothness) characteristics of the cloud model. The 
temperature of these two layers is varied by increments of 0.04 K over a range of 4 K (for a total 
of 100 experiments). We will refer to this experiment as Case 0. Figure 3.3. 1 shows the coldest 
and the warmest soundings used. All other conditions are held constant (unless stated 
otherwise), including a constant PBL height of 50 mb (-500 m) and a constant initial velocity of 
the updraft (0.2 m/s). 
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Figure 3.3.1 Skew-T diagram for environmental sounding used in the tests. The two lowest 
grid layers have a temperature that is varied between the two temperature soundings shown in 
the figure. Model levels are denoted by circles on the temperature sounding. 
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3.3.1.2.2 


Discontinuities of the Original Shallow Scheme 


To illustrate the discontinuous behavior of the original updraft scheme. Figure 3.3.2 
shows the mass fluxes in the series of 100 experiments of Case 0 at level 11 (915 mb), as a 
function of the temperature perturbation in the surface layer. Several jumps are present in the 
series, where small changes in the temperature of the initial updraft parcel produce large changes 
in the mass flux profile. Since the shallow-cloud tendencies, as well as the tendencies fed into 
the grid-resolved variables, depend directly on the mass flux profile, it appears that the current 
scheme can produce abrupt changes in the calculated grid-scale tendencies. 

Further analyses of these results have shown that the jumpiness in the mass flux is due in 
part to the large importance given to the grid layers in the scheme's computations. For example, 
the detrainment at the top of the updraft is performed over an integer number of layers (up to the 
level LTOP, which is the layer in which the cloud top is diagnosed). Therefore, when the 
number of detraining layers changes from, say, 1 to 2, the mass fluxes change abruptly and 
significantly in all layers involved. Other layers that have an important impact on the results are 
the layer of equilibrium temperature (LET) and the layer defining the thermodynamic properties 
of the starting updraft. Whenever these layers change by a unit (one whole layer), the updraft 
mass flux distribution can change drastically. 

Figure 3.3.3 shows the series of important heights diagnosed with the original updraft 
algorithm, as a function of the perturbation temperature. ZPBL is the assumed PBL height that is 
set constant at 500 m in all calculations. The fine line labeled ZLCL is the diagnosed lifting 
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Mass Fluxes at Level 1 1 



Figure 3.3.2 Diagnosed mass fluxes for the original updraft algorithm at Level 1 1 as a function 
of perturbation temperature in the initial shallow-cloud updraft parcel. UMF: mass flux at the 
top of the layer. UER: entrainment flux in the layer. UDR: detrainment flux in the layer. 
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Heights for Original Scheme 



Figure 3.3.3 Important heights diagnosed with the original shallow-cloud updraft algorithm. 
ZPBL: height of the PBL. ZLCL: lifting condensation level. CLDTOP: cloud top height. 
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condensation level that is based on the temperature and moisture of the surface layer. Here, the 
updraft release height is always at ZPBL because ZPBL < ZLCL in all of the 100 experiments of 
Case 0. Finally, the line with circles labeled as CLDTOP is the diagnosed top of the updraft (or 
cloud). The discrete behavior of the cloud top is a serious problem that complicates the 
interpretation of the scheme's results and has feedbacks to many of the shallow-cloud scheme’s 
calculations. Also, there are two regions in Figure 3.3.3 where CLDTOP drops suddenly to a 
minimum value. These regions correspond to conditions for which the algorithm diagnoses no 
shallow clouds because the condition ZLCL > CLDTOP is met. A further discussion of the 
appropriateness and implicaitons of this condition is presented below. While it is true that the 
series of LET and LTOP must jump between integer numbers (because they refer to grid levels), 
it is not necessary true that the algorithm has to reflect these discontinuities in such a strong way. 
In fact, LET and LTOP represent the integer locations of well defined heights that are continuous 
in the vertical, namely, the height of no-buoyancy (ZLET) and the height of zero vertical 
velocity for the updraft parcel (ZTOP), respectively. 

As it is now, the shallow cloud scheme is activated (i.e., produces condensed water mass, 
which is a source of NBC and defines cloud mass fluxes) depending on the relationship between 
the ZLCL and CLDTOP. If CLDTOP > ZLCL, then the cloud model is activated and mass 
fluxes are calculated. Although this condition appears reasonable, a closer look indicates that its 
consequences are not completely straightforward. Figure 3.3.4 shows the total detrained mass 
fluxes of water vapor and liquid water that are diagnosed by the original scheme. Again, the 
discreteness of the fluxes is evident in the figure. The two regions with large drops in the vapor 
detrainment fluxes in Panel a correspond to the experiments in which the shallow scheme does 
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Detrained Water Vapor Mass Flux 



Figure 3.3.4 Sources of detrained updraft water vapor and liquid water diagnosed with the 
original updraft algorithm, (a) Total detrained water vapor flux (regions with low values around 
DT = 2.5 K and DT = 4 K are not fed back to the calling model because CLDTOP<ZLCL). (b) 
Total detrained liquid water. 
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not activate because of the triggering condition discussed above (also see Figure 3.3.3). The 
detrained liquid water (that defines the source of shallow water clouds) is also discontinuous. 

Two subtle effects of the shallow-cloud triggering condition are illustrated in Figure 
3.3.4. First, 7.1 . Cl. is calculated with an undiluted parcel, but the CLDTOP calculations include 
the effects of dilution with the environment. Thus, it is possible that, although CLDTOP > 
ZLCL, the environmental dilution precludes the updraft from actually reaching saturation. In 
this condition the algorithm will produce mass flux profiles and shallow-cloud tendencies, but 
there will be no production of liquid water at all! This situation is seen in the figure in the region 
2.2 K > DT > 1.5 K, where the cloud fluxes are active, but there is no liquid-water detrainment. 

The second effect is also illustrated in Figure 3.3.4. The discrete character of CLDTOP 
produces abrupt changes in the status of the triggering condition, as can be seen in the regions 
around DT=2.2K (ON to OFF) and DT=2.55K (OFF to ON). Especially when the jump is from 
OFF to ON, the jump in the shallow-cloud parameterization tendencies can be large, because 
there already is a large updraft occurring when the ON condition is triggered. 

Another discontinuity that the updraft algorithm possesses, but which has not been 
illustrated here, is in the definition of the thermodynamic properties of the cloud-initiating 
parcel. A discrete number of model layers near the surface is considered at each time step to 
define these properties, so that as the ZPBL grows continuously, the source layer grows 
discretely. That is, the source layer is defined as two or more discrete model layers, beginning 
from the surface. Understanding and correcting this problem is more straightforward than the 
problems with the cloud model described previously. 
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In summary the problems detected in the updraft algorithm of the shallow-cloud scheme 
can lead to discontinuities in the mass flux profiles that may arise from: (a) discrete definition of 
the updraft-property layer, (b) discrete definition of the releasing layer for cloud-initiating 
parcels, and (c) discrete definitions of the LET layer and the upper moist detrainment layer. 
Furthermore, a problem with the triggering condition has been identified. The use of ZLCL in 
defining the triggering condition appears to be inappropriate and unnecessary. The improved 
updraft scheme developed in this study tries to use corresponding continuous definitions for 
these significant levels in its calculations. Results with the revised scheme are presented below 
in the next sections. 


3.3.1.2.3 Basic Results with the Improved Updraft Algorithm 


The revised shallow-cloud updraft scheme tries to retain the physical basis of the original 
Kain-Fritsch scheme and its conservative properties, while also stressing the importance that its 
diagnosed mass fluxes should be smooth and continuous. To attain this goal, the revised scheme 
requires that the environmental and/or initial updraft properties must change in a smooth and 
continuous manner. Here, we only present the results for the same experimental conditions 
(Case 0) discussed in the previous section, but using the revised shallow-cloud updraft scheme. 

Figure 3.3.5 shows the mass fluxes at level 11, using the revised updraft scheme 
(compare to Figure 3.3.2). The new results are indeed smooth and continuous. The results in the 
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Air Mass Fluxes (kg/s) 



Figure 3.3.5 Diagnosed mass fluxes using the revised updraft algorithm at Level 1 1 as a 
function of perturbation temperature in the initial shallow-cloud updraft parcel. UMF: mass flux 
at the top of the layer. UER: entrainment flux in the layer. UDR: detrainment flux in the layer. 
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figure clearly show that, as DT grows, the evolution of detrainment rate, the increase of mass 
flux, and later the development of entrainment each have a simple physical meaning that is 
consistent with the expected behavior at a fixed height when the updraft grows in depth through 
that level. 

Figure 3.3.6 shows the corresponding series of significant heights to be compared with 
Figure 3.3.3. The updraft top height, ZTOP, increases smoothly and continuously as the updraft- 
initiating parcel becomes warmer. This figure also shows some additional relevant heights that 
are diagnosed with the new scheme. Of particular interest is the height labeled as ZLCLU that 
corresponds to the height where the entraining/detraining updraft reaches saturation and 
condensed water begins to form. In this example, ZLCLU is significantly higher than ZLCL, the 
lifting condensation level based on the undiluted initial parcel (which of course is exactly the 
same as the one diagnosed by the original shallow-cloud scheme in Figure 3.3.3). Also, the level 
ZLET is the no-buoyancy level that corresponds to the LET level of the original scheme. 

Finally, Figure 3.3.7 shows the total detrained water vapor and liquid water fluxes, which 
can be compared to the corresponding fluxes of Figure 3.3.4. The smoothness of the new fluxes 
contrasts dramatically with those calculated using the original scheme. In particular, the liquid 
water fluxes are non-zero in the entire domain of this test. Furthermore, they appear to have a 
weak minimum at intermediate values of initial parcel temperature (DT - 3 K). One would 
expect that the sources and tendencies derived from the shallow convection scheme will also 
behave in a more physically plausible way. 

It is now possible to seek a physical interpretation for the behavior of the shallow 
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Height (m) 



Figure 3.3.6 Important heights diagnosed with the revised shallow-cloud updi 
ZPBL: height of the PBL. ZLCL: lifting condensation level. ZTOP is the heigl 
top. 7.1 . Cl -1 1 is the height where the updraft becomes saturated. ZLET. is the n 

level. 




Detrained Water Vapor Mass Flux - New Scheme 




Figure 3.3.7 Sources of detrained updraft water vapor and liquid water diagnosed with the 
revised updraft algorithm, (a) Total detrained water vapor flux (regions with low values around 
DT = 2.5 K and DT = 4 K are not fed back to the calling model because CLDTOP<ZLCL). (b) 
Total detrained liquid water. 
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convection scheme, which was not possible in the original scheme due to the discontinuous 
behavior of the liquid fluxes in Figure 3.3.4b. A close inspection of Figure 3.3.7b, however, still 
shows some artifacts in the experiment series of fluxes, but they are minor compared to the 
original results and should not preclude the continuation of this work. Therefore, we were able 
to move forward with a broader series of tests. 


3.3.1.2.4 Further Results with the Improved Updraft 
Algorithm 

In this section we present results from the revised shallow-cloud updraft algorithm using 
changes in other conditions affecting the cloud-initiating parcels. 


3.3.1.2.4.1 Case 1: Changes in the PBL Height 

For Case 1, the PBL height was varied from 400 to 700 m in 100 steps of 3 m each. All 
other conditions were held constant as in Case 0. (Note that the temperature perturbation of the 
surface layer in Case 1 and the following cases is fixed at 2 K, exactly in the middle of the 
sensitivity range applied in the previous section.) Changing the PBL height directly affects the 
releasing height of the updraft-initiating parcel, as well as the other properties of the initial 
parcel. The latter effect is due to the fact that the scheme assumes that the initial parcel has the 
average properties of the surface-based source layer corresponding to 20% of the PBL depth. 
Figure 3.3.8 summarizes the results for this sensitivity test for the improved (left panels) and 
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Detrained Air Mass Flux - New Scheme 



Figure 3.3.8. Results from the shallow convection as a function of change in PBL height. 
Figures 3.3.8a are computed with the revised updraft algorithm. Figures 3.3.8b are based on the 
original updraft algorithm. Upper panels are critical heights, middle panels are total detrained air 
mass fluxes, and lower panels are total liquid water fluxes. 








Figure 3.3.8b (Continued) 
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original (right panels) updraft algorithms. The upper panels show the series of heights for the 
two schemes. The old scheme maintains a constant cloud-top height until the diagnosed ZLCL 
abruptly rises at ZPBL-615 m (because an extra integer layer was considered in the surface 
source layer), causing the scheme to suddenly shut off. The revised scheme, however, diagnoses 
a gradual change in all heights, including a gradual increase of the ZLCL as the surface layer 
grows in depth along with the ZPBL. 

The middle panels show the total detrained air mass fluxes. Again, the original scheme is 
highly discontinuous and the mass flux stops abruptly. The new scheme, on the other hand, 
diagnoses a gradual decrease of the mass fluxes as the release height increases. The physical 
reason for this is that the density of the initial parcel decreases as the release height slowly rises. 
Finally, the bottom panels show the diagnosed total liquid water detraining fluxes. In this case 
the old scheme does not produce any liquid water, although the updraft scheme is active in most 
of the series. Again, this is because of the dilution effect on the parcels and the detrainment 
algorithm used by the scheme. On the other hand, the new scheme produces a continuous source 
of liquid water as long as ZLCLU is lower than ZTOP (see upper panels). When that condition is 
not true the liquid fluxes drop smoothly to zero. 


3.3.1.2.4.2 Case 2: Changes in Initial Vertical 
Speed 

In the next series of tests the initial vertical speed of the updraft parcel is varied from 0.2 
to 10.2 ms' 1 in 100 steps of 0.2 ms' 1 each, while all other conditions are held constant, as in the 
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temperature tests (Case 0). Figure 3.3.9 shows results for Case 2 based on the original and 
revised updraft algorithms. The most relevant differences produced by the improved algorithm 
are the smoothness in the liquid water fluxes and the height variations, compared to the original 

approach. 

3.3.1.2.4.3 Case 3: Changes in Mixing Ratio 
of the Surface Layer 


Finally, the tests are repeated, but now changing the water vapor mixing ratio of the 
surface source layer that helps to define the thermodynamic properties of the updraft-initiating 
parcel. All other variables in Case 3, plus the temperature of the source layer, are fixed as in the 
initial test, Case 0. The increase of the mixing ratio in the surface layer makes the LCL of the 
parcel (for both dilute and undilute versions) decrease and also induces a change in the saturation 

state of the updraft. 

Figure 3.3.10 shows the results for this series of tests. The variation of the heights is 
more complex in this case. The ZLCL, as well as the ZLCLU, decreases as the initial parcel 
becomes more humid. When ZLCL becomes lower than ZPBL, the parcels are released from the 
ZLCL. The discontinuous behavior of the liquid water detrainment mass fluxes calculated with 
the old scheme represents a serious inconsistency that is greatly reduced with the revised 
scheme. Although the improved updraft algorithm is still not completely smooth, the mass flux 
perturbations are greatly reduced so that their impacts on the NBCs and the feedbacks to the 
resolved-scale tendencies should remain negligible. 
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Detrained Liquid Water Mass Flux - New Scheme 


Initial Vertical Speed (m/s) 

Figure 3.3.9. Results from the shallow convection as a function of change in the initial vertical 
speed of the updraft parcel. Figures 3.3.9a are computed with the revised updraft algorithm. 
Figures 3.3.9b are based on the original updraft algorithm. Upper panels are critical heights, 
middle panels are total detrained air mass fluxes, and lower panels are total liquid water fluxes. 













Heights for New Scheme 



Detrained Air Mass Flux - New Scheme 




0.014 0.0145 0.015 0.0155 0.016 0.0165 0.017 0.0175 0.018 0.0185 0.019 

Water Vapor Initiating Parcel (kg/kg) 

Figure 3.3.10 Results from the shallow convection as a function of change in the mixing ratio of 
the surface layer. Figures 3.3.10a are computed with the revised updraft algorithm. Figures 
3.3.10b are based on the original updraft algorithm. Upper panels are critical heights, middle 
panels are total detrained air mass fluxes, and lower panels are total liquid water fluxes. 
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3.3. 1.3 Improvements to the Closure Algorithm 


With the improved updraft algorithm it is possible now to test the second part of the 
shallow-convection code. In this part, the closure algorithm is used to derive tendencies for the 
environmental thermodynamic profiles based on the updraft mass fluxes. The available buoyant 
energy (ABE) closure algorithm (ABE is the same as CAPE used in the Kain-Fritsch deep 
convection scheme) determines the amount of convection according to the change in ABE 
induced by the unit updraft. In the shallow convection scheme this ABE closure is supplemented 
with a TKE closure that determines the amount of convection according to the TKE in the 
boundary layer. The TKE closure is most appropriate for very shallow clouds driven primarily 
by boundary layer turbulence. For intermediate cloud depths, a hybrid regime is defined in 
which the amount of convection is interpolated between the rates calculated by the ABE and the 
TKE closures. 


3.3.1.3.1 Case 0: Results with Original Closure Algorithm 

We next applied the original closure algorithm to the same Case 0 discussed in Section 
3.3. 1.2 in which the temperature of the surface layer is varied in a set of 100 small increments of 
0.04 K. Here, we used the revised updraft algorithm from that previous section to derive the 
updraft mass fluxes that are input to the closure algorithm. Results are presented in Figures 
3.3.11 and 3.3.12. 

Figure 3.3.1 la shows the ABE calculated by the updraft algorithm for the undiluted 
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Figure 3.3.1 1 Results for Case 0 using the original closure algorithm, (a) ABE: buoyant energy 
obtained in undilute updraft algorithm. ABEG: buoyant energy calculated in the original closure 
algorithm, (b) Updraft multiplication factors calculated in the original closure algorithm. 

AINC2 is the multiplication factor calculated with the shallow TKE closure, AINC4 is calculated 
with the ABE closure, and IINC3 is the multiplication factor for the hybrid (transition) closure. 
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Figure 3 3 12 Neutrally buoyant cloud source rates, based on detriment from the shallow 
cloudupdraft calculated with the original closure algorithm for Case 0. (a) Cloud area fact, on 
source rate, (b) Liquid water source rate. 
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parcel used to compute the updraft. Also shown is the ABE diagnosed by the original closure 
algorithm (ABEG). Note that there is a large discrepancy between these two values of ABE. 

The reason is that the computation of the ABE in the closure algorithm is performed for a parcel 
starting at the undilute LCL, without regard to the fact that the parcel in the updraft algorithm 
can originate at the PBL. This condition occurs quite often, whenever the LCL is above the 
PEL. It represents a lack of consistency between the two sections of code that can have serious 
consequences because the ABE closure is highly dependent on the difference between the 
updraft ABE and the closure ABE (ABEG). In the original code both buoyant energies often are 
calculated for parcels with different properties! The apparent solution for this inconsistency is to 
use the same parcel in the closure section as it was used in the updraft scheme. Figure 3.3. 1 lb 
shows the convective amplification factors calculated by the three shallow-cloud closure 
regimes. These factors are the numbers by which the initial updraft mass flux has to be 
multiplied in order to attain closure. AINC2 is the multiplication factor calculated for the TKE 
closure, AINC4 is the factor calculated for the ABE closure, and AINC3 is the hybrid (or, 
transition) closure factor (the one that would be used in the case shown in the figure). The figure 
shows that values of AINC2 are almost constant for all 100 experiments in the senes. The values 
of AINC4 are very close to zero, but are not exactly zero. In each case the original closure 
algorithm sees an initial non-zero value of ABE coming from updraft algorithm and then 
calculates a new ABEG which is zero (because a different parcel was used). In effect, the 
closure algorithm thinks that it has overestimated the convective factor, AINC4, and therefore 
reduces the factor by 50%. The procedure is repeated until the maximum number of iterations is 
reached (10) and that value is selected as the final AINC. The transition factor for the hybnd 
closure, AINC3, is close to AINC2 with a decreasing tendency as the temperature perturbation 
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increases. 


Figure 3.3.12 shows the source rates for the key NBC (neutrally buoyant cloud) variables 
diagnosed by the original closure scheme. Figure 3.3.12a shows the source rate of cloudy area 
fraction and Figure 3.3.12b is the source rate for liquid water. The series shows the rates for 
each of the different layers at which updraft mass is detrained, but the individualization of each 
layer is irrelevant for the present discussion. The salient feature of Figure 3.3.12a is the 
discontinuous series of the cloudy area source term. With all the effort put into producing 
continuous updraft fluxes, one wonders how these source rates can collapse abruptly to zero. 

The reason lies in the algorithm that computes the tendencies of the NBC variables. We saw in 
the section discussing the updraft algorithm that an updraft can exist without producing any 
liquid water. That means that one can have a non-zero detrainment mass flux with a zero liquid- 
water detrainment rate. The original code assumed that any detrainment flux produces a cloudy 
area fraction, which is not the case when the detrainment has no liquid water. As the liquid 
water source term in Figure 3.3.12b goes to zero the cloud-area source term may not go to zero, 
which means that a "dry" detrainment is being computed as a source for the cloudy area! The 
solution to this problem requires just to keep track of the effective air mass flux that is detrained 
when liquid water is present and to use only this portion to compute the cloud area source term. 

3.3.1.3.2 Case 0: Results with Revised Closure Algorithm 

Figures 3.3.13 and 3.3.14 show results obtained with the revised closure algorithm 
described above. The ABE series shown in Figure 3.3.13a is the same as that shown for the 
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Figure 3.3.14 Neutrally buoyant cloud source rates, based on detrainment from the shallow 
cloud updraft calculated with the revised closure algorithm for Case 0. (a) Cloud area fraction 
source rate, (b) Liquid water source rate. 
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original scheme in Figure 3.3.1 la, because both are based on the same shallow-cloud updraft 
scheme. However, although the final ABEG series is lower than ABE, it is still non-zero. The 
similarity in the shapes of the ABE and ABEG series in Figure 3.3.13a is an indication that the 
same parcel properties and updraft algorithm are being used to calculate ABE in the updraft 
algorithm and in the closure algorithm. Meanwhile, Figure 3.3.13b shows the convective factors 
calculated with the revised closure scheme. The TKE convective multiplication factor, AINC2, 
is similar to that obtained with the original scheme. However, the ABE convective factor, 

AINC4, is substantially different. Instead of remaining zero it increases rapidly as the 
temperature perturbation increases, from zero to the maximum allowed value of 1000. Thus, the 
resulting hybrid-closure multiplication factor, AINC3, increases slightly instead of decreasing, as 
occurred with the original scheme (compare to Figure 3.3.1 lb). The fact that the convective 
fluxes grow when calculated according to the revised closure algorithm as the initial updraft 
parcel becomes warmer has more physical appeal than the opposite case found with the original 

closure algorithm. 

Next, Figure 3.3.14 shows the source rates for the NBC determined from the shallow- 
cloud updraft variables, as calculated according to the revised closure algorithm. Now the 
cloudy-area tendencies in Figure 3.3.14a are continuous and go to zero as the liquid-water source 
rates go to zero in Figure 3.3.14b. Again, this is physically more consistent behavior than the 

results obtained with the original algorithm. 

At this point it is illustrative to check the updraft and closure algorithms to see exactly 
what type of thermodynamic profiles they produce and to see how the reduction in ABE is being 
accomplished. Figure 3.3.15 shows a Skew-T diagram with results based for Case 0 and for the 
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Figure 3.3.15 Skew-T diagram for results with revised algorithms in Case 0, shown for the 
largest surface temperature perturbation (4 K). A description of the lines is given in the text. 
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maximum perturbation in the surface layer temperature (i.e., DT=4 K). The bold lines in the 
figure represent the environmental temperature and dew point profiles. The fine line soundings 
correspond to the final soundings calculated with the closure scheme (i.e., the soundings towards 
which the closure scheme will calculate its tendencies). The short fine lines correspond to the 
parcel trajectories in temperature and dew point for the shallow-cloud updraft parcel calculated 
with the updraft algorithm. Several features are to be noted in this figure. First, the ABE is 
calculated (approximately) as the positive area between the parcel temperature path and the 
environmental temperature sounding. We see then that the final environmental sounding indeed 
has a lower ABE than the initial sounding. The decrease in ABE, however, is not large. Also, 
just by looking at the initial and final environmental temperature profiles, one can see an upper 
level cooling and a lower level warming that are induced by the convective scheme. This 
suggests a tendency for the profile to destabilize due to the convection! In other words, by 
reducing the ABE relative to the updraft parcel, the instability of the environmental sounding 
actually is increased. This is a very interesting result and appears to match the expected behavior 
of shallow cloud fields that, in certain circumstances, they can tend to destabilize the initial 
environment. Given that the initial temperature profile in these experiments is isothermal above 
the boundary layer (i.e., very stable) and that convection should move the thermal profile 
towards a moist adiabatic condition, the result in Figure 3.3.15 is indeed consistent with 
expectations. 

A final feature to observe in this diagram is the tendency for the environmental soundings 
below the shallow-cloud updraft. Strong drying and cooling can be observed in the region just 
below the updraft base near 960 mb. The reason for this is that the source of the updraft mass 
flux is located there, but the thermodynamic properties of the shallow-cloud updraft are taken 
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from the surface layer (985-1000 mb). Therefore, when the surface layer is moister and warmer 
than the mass-source layers, the energy and water budget equations can induce strong tendencies 
in the mass-source layers. This behavior may not be physically realistic. Perhaps a better 
approach would be to assign the full boundary layer (BL) to be the source of the updraft mass 
flux. In this way the effect on the tendencies would be reduced and probably would be smoothed 
out quite easily later on by the mixing forced by the model's turbulence parameterization. 


3.3.1.3.3 Further Results with Revised Closure Algorithm 
3.3.I.3.3.I. Revised Closure Algorithm in Case 1 

We present here results of an experiment series based on the revised shallow-cloud 
updraft and closure algorithms for Case 1, in which the PBL height was increased in steps of 3 m 
from 400 to 700 m. Figure 3.3.16 shows the series of ABE (Fig. 3.3.16a), convective 
multiplication factors (Fig. 3.3.16b), and the NBC source rates for cloud area and liquid water 
content (Figs. 3.3.16c and d). For lower parcel releasing heights (i.e., at the lower values of the 
BL height), the parcels have non-zero ABE that decrease slowly as the PBL grows. Figure 
3.3.16b shows that the convection multiplication factors diagnosed by the revised closure 
algorithm also decrease slightly for shallow PBLs. For larger PBL heights (>550 m) the initial 
ABE drops rapidly to zero and the convective factors also decrease, especially AINC4, which 
comes from the ABE closure. The hybrid-closure multiplication factor, AINC3, is only slightly 
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Figure 3.3.16 Results for the shallow convection scheme with the revised closure algorithm in 
Case 1. (a) Initial and final buoyant energies, (b) Convective multiplication factors. 
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Figure 3.3.16 (Continued) Results for the shallow convection scheme with the revised closure 
algorithm in Case 1. (c) Cloudy area source rates, (d) Liquid water source rates. 
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reduced because it is dominated by the TKE closure, which is more or less constant in this case. 

Figure 3.3. 17 shows Skew-T diagrams for the lowest and highest PBL heights considered 
in the Case 1 computations. In the first case (left panel) the initial ABE is related to the positive 
area between the updraft path and the original temperature profile. This ABE decreases in the 
target environmental profiles as there is warming in the lower levels of the updraft layer. The 
updraft parcel path can be seen to reach saturation at -900 mb. Intense drying and cooling is 
seen at -960 mb (below the updraft base), which once again is related to the definition of the 
updraft mass-source layer. In the right panel of Figure 3.3.17, results for the highest PBL height 
are illustrated., Here, the updraft parcel is everywhere cooler than the environment, which 
explains the zero ABE diagnosed in this case (see Fig. 3.3.16a). The TKE closure still diagnoses 
a convective flux that slightly modifies the environment, even though the updraft fails to reach 
saturation. 


3.3.1.3.3.2 Revised Closure Algorithm in Case 2 


Results for Case 2, in which the initial vertical velocity of the shallow-cloud updraft 
parcel is incremented between 0.2 and 10.2 ms , can be seen in Figures 3.3.18 and 3.3.19. The 
series of ABE in Figure 3.3.18 shows that the final ABEG is always lower that the initial ABE, 
and in most cases it is zero. The convective multiplication factors in Figure 3.3.18b show that 
AINC4 decreases as the initial vertical velocity increases. This is natural because, as the initial 
speed of the parcel increases, so also does the mass flux of the "unit updraft. Thus, the final 
updraft mass flux corresponds to a decreasing amplification factor for the initial updraft. 
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Figure 3.3.17 Skew-T diagram for shallow convection results with the revised 

algorithms for Case 1. (a) Using lowest PBL height, (b) Using highest PBL height. 
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Figure 3.3.17 (Continued) Skew-T diagram for shallow convection results with the 

revised algorithms for Case 1. (a) Using lowest PBL height, (b) Using highest PBL height. 
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Figure 3.3.18 
algorithm in Case 2. 


Results for the shallow convection scheme with the revised closure 
(a) Initial and final buoyant energies, (b) Convective multiplication factors. 
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Figure 3.3. 18 (Continued) Results for the shallow convection scheme with the revised closure 
algorithm in Case 2. (c) Cloudy area source rates, (d) Liquid water source rates. 
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Figure 3.3.19 Skew-T diagrams for results of shallow convection scheme with the revised 
algorithms for Case 2. (a) Using smallest initial vertical speed, (b) Using maximum initial 
vertical speed. 
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Figure 3.3. 19 (Continued) Skew-T diagrams for results of shallow convection scheme with 
the revised algorithms for Case 2. (a) Using smallest initial vertical speed, (b) Using maximum 
initial vertical speed. 
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Figure 3.3.19 shows the Skew-T diagrams corresponding to the lowest (0.2 ms' 1 ) and 
highest (10.2 ms' 1 ) initial vertical velocities. In Fig. 3.3.19a we see that there exists some initial 
ABE that decreases as the environment is modified by convection. The updraft reaches 
saturation in the upper levels near 900 mb. In Figure 3.3.19b for the highest vertical velocities, 
the induced tendencies are much greater than for the other case with low vertical velocities. The 
final ABE is zero (when the virtual effect of water vapor is considered). The updraft extends 
through a deep layer (up to -860 mb), where it is saturated and where liquid water is being 
detrained. The drying of the updraft mass-source layer below the updraft base in this case is very 
large and most probably is unphysical. 

3.3.1.3.3.3 Revised Closure Algorithm in Case 3 

Results obtained with the revised updraft and closure algorithms for Case 3 are shown in 
Figures 3.3.20 and 3.3.21. Figure 3.3.20a shows that as the water vapor content of the shallow- 
cloud updraft parcel increases, the convective multiplication factor for the ABE closure increases 
rapidly. Due to the TKE dominance in the hybrid (transition) closure, however, the ABE is only 
partially reduced in the final environmental profile (Fig. 3.3.20b). As the initial parcel becomes 
saturated, the induced tendencies grow substantially. This is especially evident in the Skew-Ts 
presented in Figure 3.3.21. Note that in the right panel of that figure, the parcel trajectory is 
saturated aloft from the shallow-cloud updraft base (where the temperature and dew point lines 
overlap at -965 mb). The temperature sounding of the parcel trajectory is, however, 
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Figure 3.3.20 Results for the shallow convection scheme with the revised closure algorithm in 
Case 3. (a) Initial and final buoyant energies, (b) Convective multiplication factors. 
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Figure 3.3.20 (Continued) Results for the shallow convection scheme with the revised closure 
algorithm in Case 3. (c) Cloudy area source rates, (d) Liquid water source rates. 
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Figure 3.3.21 Skew-T diagram for results with revised algorithms in Case 3. (a) Using the 
smallest water vapor surface perturbation. A description of the lines is given in the text 
following Figure 3.3.15. 
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significantly different than the moist adiabat. This effect is due to the entrainment of 
environmental air into the cloud updraft. A test of this case with entrainment set to zero yields, 
quite naturally, a moist adiabat trajectory on the Skew-T (not shown). 

3.3.1. 3.3.4 Extension of the Mass-Source Layer to 

the Full PBL 

As pointed out previously restricting the shallow-cloud updraft mass-source layer to the 
upper levels of the PBL, while taking its thermodynamic properties from near the surface, can 
induce large tendencies in the region of the mass source layer. These large tendencies appear to 
be unphysical. Therefore, a test was conducted by modifying this algorithm to redefine the 
mass-source layer as the complete PBL. Figures 3.3.22 and 3.3.23 show results for Case 0 with 
this modification of the source-layer algorithm. In the shallow-cloud updraft layer above 950 mb 
(in Figure 3.3.23), the results are similar to those obtained previously (in Figure 3.3.15), but the 
tendencies in the PBL shown in the Skew-T diagram of Figure 3.3.23 are much smaller than in 
the previous Case 3 experiment in which the mass-source layer was restricted to the upper 
portion of the PBL. 
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Figure 3.3.22b Results for the shallow convection scheme with the revised closure 

algorithm in Case 0 with the mass-source layer extended to the entire PBL. (c) Cloudy area 
source rates, (d) Liquid water source rates. 
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Figure 3.3.23 Skew-T diagram for results of the shallow convection scheme with revised 
algorithms in Case 0, but for the mass-source layer extended to the full PBL. The results are 
based shown for the largest surface temperature perturbation (4 K). A description of the lines is 
given in the text following Figure 3.3.15. 
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3.3.2 Calibration and Sensitivity Tests 


3.3.2.1 Introductory Remarks 


We describe here improvements and sensitivities of the Penn State (PSU) shallow 
convection parameterization. The original parameterization was developed in Deng (1999) and 
has been described in Deng et al. (2002a, 2002b). In Munoz et al. (2001) the scheme was 
applied in the 3-D MM 5 model Version 3.3 to simulate a summertime convection case over the 
central U.S. (see Section 4.3). Although the results were encouraging because the shallow 
convection scheme was able to partially capture the spatial distribution of shallow clouds, they 
also showed some systematic underprediction of cloudiness and a general patchiness in the cloud 
fields, that prompted the revisions and improvements described here. 

The PSU shallow convection scheme can be considered an extension of the Kain-Fritsch 
(KF) deep convection parameterization scheme (Kain and Fritsch, 1990) to account for the 
effects of shallow non-precipitating clouds. As a standard convection parameterization scheme, 
the PSU shallow convection scheme consists of a triggering condition, a cloud model and a 
closure assumption. Because it is designed to represent shallow convection, however, the new 
scheme differs from other convection parameterizations in several respects. On the one hand 
some important elements of deep convection parameterizations lose their importance when 
applied to shallow convection, such as the triggering condition and the downdraft sub-model. 
The first has less importance because shallow clouds are much more numerous in a given area 
than deep-convection events. Also, since shallow clouds are by definition non-precipitating (or 
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nearly so), they do not require consideration of downdrafts that are formed by the evaporation of 
precipitating hydrometeors. 

On the other hand, the representation of shallow clouds in mesoscale models poses some 
problems that are not as relevant in the deep convection case. For example, it is important to call 
the parameterization scheme much more often than the deep-convection routines, since the 
shallow clouds interact rapidly with many of the other physics of the calling model, mainly the 
land-surface, turbulence and radiative parameterizations. If the shallow convection is called 
infrequently, these interactions may be lost or distorted. Another important special feature of 
shallow clouds is that they have much smaller thermodynamic perturbations, compared with the 
deep-convection case. This characteristic makes the definition of the initial properties of the 
updraft parcel a very sensitive part of the parameterization. It also demands that updraft and 
mixing algorithms that have proven useful in the parameterization of deep convection be re- 
evaluated in the more subtle conditions occurring in shallow convection, where the buoyancy or 
liquid water content of the updraft may differ only marginally with respect to the environment. 

A final special feature of the shallow-convection scheme is its treatment of the liquid 
water generated in the updraft. Unlike the Kain-Fritsch deep-convection scheme, where 
detrained liquid water is fed back immediately into the resolved water variables, the shallow 
scheme retains the detrained condensed water in sub-grid cloud variables (area and water content 
of sub-grid neutrally buoyant clouds, NBCs). The water mass of the NBCs is gradually put back 
into the resolved-scale variables by means of mixing and other parameterized processes. 

This section describes further tests and improvements of the shallow-convection scheme, 
in addition to those discussed in Section 3.3.1, aimed at making it more robust and useful in the 
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parameterization of shallow convection for mesoscale and climate models. It also validates the 
improved scheme against a published data set from a trade-wind shallow convection case. 

3.3.2.2 New Features in the Shallow Convection Parameterization 
3.3.2.2.1 Parcel Properties 

The original shallow scheme defines the thermodynamic properties of the updraft 
initiating parcel based on the average properties of the lowest 20% of the planetary boundary 
layer (PBL). The starting height of the updraft is defined as the minimum between the PBL 
height and the lifting condensation level (LCL) of the initial updraft parcel. The results shown in 
this section have been obtained with a simpler definition of the initial parcel properties. The 
buoyancy at the top of the PBL is calculated and a given perturbation to this value is assigned to 
the parcel, which is also considered to be saturated, but with no liquid water content. The 
buoyancy perturbation is a parameter of the scheme, as well as the perturbation of the initial 
vertical velocity. As we show later, very small buoyancy perturbations are sufficient to create a 
shallow-cloud layer. This new scheme is simpler than the original scheme and concentrates on 
the important physical variable (initial buoyancy) more than on a prescribed fraction of the PBL. 
Note that this new approach does not imply that the PBL top and cloud layer base are forced to 
coincide, since entrainment of dry air can preclude the parcel from producing liquid water near 
the updraft’s base. 


3.3.2.2.2 Updraft Model 
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Two problems have been found when applying the Kain-Fritsch cloud model in the 
shallow-convection scheme. These problems are not significant for the use of the K-F cloud 
model in the case of deep-convection parameterization, but they become important when used in 
the context of shallow convection, as was explained in the introduction. 

The continuous use of the shallow-convection parameterization (i.e., it being called every 
time step) uncovered the degree of discretization inherent in the original K-F updraft scheme (see 
Section 3.3.1). The updraft is discretized on the same vertical grid layers as in calling model, 
which can have serious consequences. For example, in the time step when the detrainment layer 
at the top of the updraft includes an additional layer, the updraft mass fluxes can decrease 
abruptly by a significant percentage (from 100% to 50% in the case that the detrainment layer 
comprises initially one and then two model layers). These sudden jumps in the mass-flux 
profiles produce noise in the time evolution of the shallow-convection tendencies, which makes 
interpretation of the results more difficult and precludes attainment of quasi-equilibrium 
conditions. 

The proposed solution to this discretization problem is to make the updraft algorithm able 
to recognize when a given model layer includes a significant transition level. We define four 
types of significant transition levels where 

1. buoyancy changes sign, 

2. liquid water becomes positive, 

3. vertical velocity becomes zero. 

4 mass flux becomes zero. 
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Once the algorithm recognizes the existence of a significant transition, it iterates inside the 
model layer to find its exact location, and then inserts this height into the updraft vertical layers. 
The transitions of types 3 and 4 also define the top of the shallow-could updraft. The logic of the 
new algorithm is somewhat more complex than that of the original scheme, but it results in an 
updraft that evolves smoothly as the environment and the initial parcel properties gradually 
change. 

The second problem relates to the small thermodynamic perturbations associated with the 
shallow clouds. Figure 3.3.24 is taken from Kain and Fritsch (1990), and it shows the conceptual 
scheme used to define the detrainment and entrainment fractions in the original K-F updraft 
scheme. The updraft parcel is mixed in different proportions with an environmental parcel. In 
the figure the buoyancy of the resulting mixed parcel is plotted against the fraction of 
environmental air used. The point of crossover, where the mixed parcel has zero buoyancy with 
respect to the environment, is the critical mixing fraction that is used later to define the 
entrainment and detrainment rates. The larger the critical mixing fraction is, the larger the 
entrainment, and vice-versa. 

This conceptual scheme is appealing and has been applied for a long time in the K-F 
deep-convection parameterization. Note however how the virtual temperature perturbation of the 
updraft in Figure 3.3.24 is rather large, as well as its water content that controls the maximum 
virtual cooling of the mixed parcels. In shallow clouds the thermodynamic perturbations are 
smaller, as well as the liquid water content, and this mixing algorithm tends to produce 
entrainment and detrainment that are mutually exclusive (i.e., all entrainment or all detrainment 
in a given layer). In the limit, for example, when the updraft parcel has no liquid water content, 
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Figure 3.3.24. Virtual temperature of a parcel resulting from the mixture of updraft and 
environmental air (taken from Kain and Fritsch, 1990). 
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the scheme has a flip-flop behavior where full entrainment is diagnosed if the updraft is buoyant 
or full detrainment is diagnosed if it is not. Conceptually the applicability of the original K-F 
mixing algorithm near cloud base is questionable, since there has been little time for mixing to 
take place, so that the buoyancy of mixed parcels may not be very relevant to the definition of 
the entrainment and detrainment profiles. Along these lines, the change to the algorithm that we 
test here is to alter the definition of the critical mixing fraction in regions near cloud base (which 
mainly affects shallow clouds). In these regions the critical mixing fraction is defined solely on 
the basis of the buoyancy of the updraft parcel with respect to the environment, using a smooth 
function of the form: 


V F_ = 


T -T 

v,u vji 

2DTM 


-,0<4' <1 

7 m 


(3.4) 


where T v>min = T v ,env - DTM, T v>u is the virtual temperature of the updraft parcel and T v , en v is 
the virtual temperature of the environment. 

The actual mixing fraction used in the entrainment/detrainment diagnosis is a weighted 
average of the K-F mixing fraction and the new critical mixing fractions, in the form: 


Z z 

xp = g'DELZSXp^ + (1 _ e'DELZS ^ 


(3.5) 
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where 'F kf is the critical mixing fraction diagnosed with the standard KF algorithm. The two 
new parameters introduced by this scheme are DTM and DELZ, and are subject to sensitivity 
evaluations in the next section. 

3.3.2.2.3 Closure 

In the results presented here we restrict the scheme to use only the TKE closure algorithm 
described in Deng et al. (2002a, 2002b), with no use of the transitional closure regime between 
TKE and CAPE closures. This is so because we intend to test the scheme in a purely shallow- 
convective situation. Evaluation of the transitional closure algorithm will be performed in future 
work. 


3.3.2.3 A trade wind boundary layer case 

As a reference for the calibration and sensitivity analysis presented later we use a trade 
wind boundary layer case described by Siebesma and Cuijpers (1995). They showed 
observations and LES results for an undisturbed period of the BOMEX experiment. Figure 
3.3.25 shows vertical profiles of mass flux (Panel a), entrainment rate (Panel b) and detrainment 
rate (Panel c). The LES results correspond to averages for the last 4 h of a 7-h simulation, and 
they report results obtained with two different methods of defining the updrafts. The case is 
forced by prescribed subsidence and radiative-cooling profiles, constant humidity advection, and 
constant sea-surface temperature. In our model calculations we use the same forcings and the 
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Figure 3.3.25. Mass flux profiles for BOMEX trade wind case, a) Mass flux profile, 
b) Entrainment rate profile, c) detrainment flux profile. Bold lines are profiles derived from 
budgets of observed data, and other two lines are derived from LES results with two methods of 
defining updrafts, (taken from Siebesma and Cuijpers, 1995). 
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same initial soundings. One of the main conclusions in Siebesma and Cuijpers (1995) is that, 
compared to LES results, current mass-flux parameterization schemes appear to underestimate 
the mixing between updrafts and environment by an order of magnitude. 

3.3.2.4 Parameter Sensitivities 

The parameter sensitivities described in this section are based on an experiment called 
Case 08, the parameters for which are shown in Table 1. 


Table 1. Parameter values used in the basic sensitivity cases Case 08 and Case 10. 


Parameter 

Case08 

Case 10 

DTVI 

0.10 


ENTPAR 

0.15 


CLOSVAL 

20 

20 

DELZS 

300 

500 

DTM 

0.20 

0.25 

D1FFK 

0.00001 

0.00005 


Figure 3.3.26 shows sensitivities of the new scheme to DTVI, the initial buoyancy 
perturbation of the updraft parcel (panels a, b, and c), and to ENTPAR, the entrainment 
parameter defined by the K-F standard scheme (panels d, e, and f). The vertical profiles shown 
correspond to 4-hour averages for mass fluxes (left panels), entrainment rates (center panels), 
and detrainment rates (right panels). 
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Figure 3.3.26. Sensitivity of new shallow convection scheme to initial parcel buoyancy 
perturbation DTVI (upper panels) and KF entrainment parameter ENTPAR (lower panels), a and 
d) are mass flux profiles, b and e are entrainment rate profiles, c and 0 are detrainment rate 
profiles. 


131 




The sensitivities for DTVI include initial parcels that are negatively buoyant (DTVI=-0.2 
K), neutrally buoyant (DTVI=0.0 K), and positively buoyant (DTVI=+0.2 K). The negatively 
buoyant initial parcels produce a shallower updraft layer and weak mass flux profiles dominated 
entirely by detrainment. The other two cases show similar mass flux profiles, except for larger 
entrainment rates near the updraft base in the case of the buoyant initial parcel. As the initial 
buoyancy increases, the mass flux near the base increases and the base height itself decreases, 
possibly due to the effect of the induced subsidence on the boundary layer growth. 

Sensitivity of the new scheme to the K-F entrainment parameter is shown in the bottom 
three panels of Figure 3.3.26. Here, as ENTPAR increases from 0.03 m/Pa (the default value in 
K-F) to 0.20 m/Pa the entrainment and detrainment fluxes increase considerably. The case with 
lower entrainment also produces a deeper updraft. 

The three upper panels of Figure 3.3.27 show the sensitivity of flux profiles for changes 
in the closure parameter, CLOSVAL. This parameter controls the magnitude of the mass flux at 
the updraft base when the TKE closure assumption is used. Panel a) clearly shows the influence 
of this parameter on the base mass flux. Increased cloud-base mass flux produces entrainment 
and detrainment profiles consistently larger, although the updraft top is not very sensitive to this 
parameter. The updraft base, however, decreases as the mass flux increases. 

The three lower panels in Figure 3.3.27 show the sensitivity of the scheme to DELZ, the 
scale height that controls the relative importance of the original K-F mixing fraction and the new 
definition of the critical mixing fraction (see section 3.3.2.2.2). This parameter greatly 
influences the entrainment/detrainment ratio in the middle levels of the updraft, where the new 
buoyancy-based critical mixing fraction tends to increase the magnitude of the entrainment and 
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Figure 3.3.27. Sensitivity of new shallow convection scheme to closure parameter CLOSVAL 
(upper panels) and height scale DELZS (variables in lower panels are defined as in Figure 
3.3.26). 
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decrease the detrainment. As a result the updraft layer is deeper and the detrainment increases in 
the updraft top layers. 

Figure 3.3.28 shows sensitivities of the new scheme to the parameter DTM. This 
parameter controls the smoothness of the buoyancy transition. As DTM decreases the transition 
becomes steeper and the critical mixing fraction jumps faster from 0.0 to 1.0 (see Section 
3.3.2.2.2). Panel a) in the figure now shows the average cloud-area fraction calculated for three 
values of DTM. Panels b) and c) show mean entrainment and detrainment rates profiles, 
respectively. Panels d), e), and f) show the standard deviation of the corresponding profiles for 
the same 4-h period in which the averages where performed. It is hard to find a systematic 
sensitivity in the mean profiles, but the influence of this parameter in the standard deviation is 
noticeable. In fact, as DTM becomes smaller the variability of the mass flux profiles increases 
significantly. This variability represents the abrupt changes from which the diagnosed updraft 
suffers on a time-step-to-time-step basis, as the buoyancy of the updraft parcel changes slightly. 
Note that the variability in the cloud-area fraction profile (note shown) does not change 
significantly, because it can be considered as a time integrated variable of the NBC, rather than 
an instantaneous function of the updraft characteristics. 

3.3.2.5 Validation and Physical Sensitivities 

Based on the sensitivities shown previously, a suitable set of parameters was selected that 
produces mass flux profiles and cloud area fractions that resemble those presented by Siebesma 
and Cuijpers (1995). This set of parameters is shown in Table 1 in the column labeled Case 10. 
For this case we performed sensitivity analyses with respect to changes in the prescribed 
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forcings. We report here sensitivities due to the specification of the subsidence profile and the 
sea-surface temperature. 

The three upper panels of Figure 3.3.29 show changes in the vertical profiles of cloud 
area, entrainment rates and detrainment rates when the subsidence profile of the basic case (Case 
8) is multiplied by a constant factor of 0.5, 1.0, and 2.0. For increased subsidence the updraft top 
is lower and the updraft layer becomes shallower. Cloud cover and mass fluxes close to the 
updraft top are greatly affected by changes in the subsidence, but are not affected so much near 
updraft base. 

The three lower panels of the same figure show sensitivities to sea-surface temperatures 
that vary -1, 0 or +1 K with respect to the basic case (Case 8). Larger sea-surface temperatures 
increase cloud area fractions, updraft depth and entrainment and detrainment fluxes. 

We finally illustrate in Figure 3.3.30 the effects of the shallow-convection scheme on the 
thermodynamic sounding of the environmental column. Panel a) is the initial Skew-T sounding 
for Case 10. Panel b) shows the sounding after seven hours when the shallow-convection 
scheme was not used, and Panel c) is the corresponding 7-h sounding when the shallow 
convection scheme was active. Among the several effects of the shallow-convection scheme on 
the thermodynamic profile we can mention: (1) shallow convection precludes the top of the 
boundary layer of becoming saturated, (2) the cloud layer becomes more humid, especially in its 
upper part, (3) the trade-wind inversion at the top of the cloud layer is more highly intensified 
(larger temperature gradient) when using the shallow-convection scheme, and 4) the top of the 
mixed layer develops a small temperature inversion. All of these effects appear physically 
realistic and consistent with the conceptual understanding of the trade-wind boundary layer. In 
the case that did not use the shallow convection parameterization, the top of the PBL becomes 
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Figure 3.3.29. Sensitivity of new shallow-convection scheme to forcings. Upper panels show 
sensitivities to subsidence profile, and lower panels are sensitivities to sea surface temperature. . 
Panels a) and d) are cloud area profiles, Panels b) and e) are entrainment rate profiles, and Panels 
c) and f) are detrainment rate profiles. 
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Figure 3.3.30. a) Initial sounding for trade-wind boundary layer case, b) Sounding after 7-h run 
without using the shallow convectionl parameterization, c) Same as b), but with the use of the 
shallow-convection parameterization. 
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fully saturated (stratus development), even though the 1.5-order TKE-predicting turbulence 
scheme has been corrected for the effects of saturation on the buoyancy production of TKE. 

3.3.2.6 Conclusions 

We have described here several improvements to the PSU shallow-convection scheme 
that make it more robust for its use in mesoscale models. The changes are demanded by the 
special features of shallow convection, especially because of the more subtle thermodynamic 
perturbations associated with shallow clouds that require more frequent use of the 
parameterization. The improvements have been made mainly in the Kain-Fritsch updraft 
algorithm that also is used in the shallow-convection parameterization. 

Use of the improved scheme in the modeling of a trade-wind boundary layer case shows 
that the entrainment parameter of the K-F updraft model may be significantly smaller than 
required to reproduce the characteristics of observed and LES-derived mass flux profiles, as 
presented in Siebesma and Cuijpers (1995). This conclusion is consistent with that reported by 
the authors of the same paper. Moreover, with the correction of the entrainment parameter, the 
new shallow-convection scheme appears to perform reasonably well for reproducing the mass 
flux profiles, while showing physically realistic sensitivities to its important parameters and to 
changes in the prescribed forcings. 

A final point of discussion concerns the shape of the detrainment profile obtained with 
the shallow-convection parameterization. Generally, it proved difficult to produce detrainment 
profiles that had maximum values close to cloud base and decreased with height, as it is typically 
the case in observed or LES-simulated shallow-cloud layers. In most of the results presented 
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here the detrainment profiles have a tendency to be uniform or increase with height. That is, the 
detrainment profile tends to peak in the aptly named detrainment layer at the top of the updraft. 
This appears to be a consistent shortcoming of many one-cloud convection schemes. However, 
it is interesting to note that in Panel c of Figure 3.3.28, the detrainment profile obtained with 
DTM=0.05 K has a shape that indeed tends to decrease with height. As discussed in Section 
3 3.2.4, a small value of DTM corresponds to a more intermittent updraft. This observation 
immediately suggests that a possible way of improving the performance of a one-cloud 
convection scheme like this may be to add an additional degree of variability to the properties of 
the initial updraft parcel on successive time steps by changing its buoyancy, vertical velocity or 
mass flux according to a prescribed distribution function. The diagnosed updrafts on successive 
time steps would then show "noise", but this will be controlled by known and physically 
consistent causes, more than by numerical or algorithmic reasons. Over a number of time steps, 
the resulting forcing of the environment by the shallow convection would have the effect of an 
ensemble of shallow-cloud sizes without having to calculate the effect of the entire ensemble at 
every time step. This approach could be computationally simple and efficient, and hence it may 
be beneficial for the representation of the overall process of shallow convection. Further 
investigations along this line will be pursued in the future. 
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4. EXPERIMENTATION IN THE 3-D MM5 


4.1 TKE-Predicting Turbulence Scheme in the 3-D MM5 

4.1.1. Introductory Remarks 

One area that is often overlooked in numerical weather prediction is proper 
treatment of turbulence in a saturated environment. The Blackadar planetary boundary layer 
(PBL) scheme (Grell et al. 1994, Zhang and Anthes 1982), for example, often produces dry- 
adiabatic or other absolutely unstable profiles in saturated layers. Although moist corrections for 
the vertical diffusion in cloudy layers have been added to the Richardson number (which appears 
in the computation of the mixing coefficient within the scheme), these problems have persisted 
and still occur quite regularly. This section presents some important considerations for 
representing atmospheric turbulence in layers saturated with respect to liquid or ice. 

4.1.2 The TKE-Predicting Sub-Model for Turbulence 

The 1.5-order turbulence sub-model that has been developed at Penn State (Gayno 1994, 
Gayno et al. 1994, Shafran et al. 2000) is now an option in version 3 of the MM5. The original 
version of this scheme, also known as the Gayno-Seaman scheme, is based on the work of 
Musson-Genon (1987) and Ballard et al. (1991). Although the scheme also contains a subgrid 
fog parameterization, it will not be discussed here and it is not included in the V3 release. The 
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fog module is still under development and will be included in a future release of the PSU TKE 
sub-model. 

Although TKE-based schemes are often called PBL schemes, they are actually turbulence 
parameterizations for the entire atmosphere. The so-called “1.5-order closure” (Mellor and 
Yamada 1974, Yamada 1977) includes a prognostic second-order equation for only the 
momentum variance, or the turbulent kinetic energy (TKE). It has been shown that 1.5-order 
schemes are reasonably accurate and more computationally efficient than the complete second 
order schemes. The remainder of this section presents some special considerations for 
representing turbulent and microphysical processes in saturated layers, including important 
improvements to the scheme beyond that described in Shafran et al. (2000). 


4.1.2.1 The TKE Equation 

As described in Shafran et al. (2000), a prognostic equation for TKE, comprised of 
source and sink terms due to buoyancy, shear, transport and dissipation, is integrated forward in 
time, along with the other mixing variables. (We have recently included the advection of TKE, 
which may be important for accurate simulation of the PBL on fine grids in transition zones 
between land and water, dry land and moist land, etc.) The vertical mixing coefficients are 
proportional to the predicted TKE profiles in time and space. The mixed-layer height is also 
diagnosed from the TKE profiles using minimum TKE threshold conditions. 

4.I.2.2. The Conservative Mixing Variables 
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At the beginning of an advection time step, the surface fluxes are computed using the 
same Monin-Obukhov similarity parameterization (Slab) used in the Blackadar scheme and 
described by Zhang and Anthes (1982). The MM5 temperature and moisture fields are then 
transformed into variables that are conservative during condensation and evaporation processes. 
The TKE sub-model converts the potential temperature (9 ) , vapor mixing ratio Oyjand liquid 
water mixing ratio (< 7 ,) into liquid water potential temperature ( 0 ,)and total water mixing ratio 
(< 7 ,) , following Betts (1973) and given by 

=*-“*■*/ (4.1) 

T c p 

<7, = <7v + . (4.2) 

where T is temperature, L Iv is the latent heat of vaporization, and e^is the specific heat at 

constant pressure. The mixing is then performed on the liquid water potential temperature and 
total water mixing ratio, and these fields are converted back to the native MM5 variables at the 
end of the time step 

Although the vertical mixing of conservative variables helps to reduce the occurrence of 
unrealistic thermal profiles in saturated layers, moist- adiabatic cloud layers are still not assured. 
We realized that the buoyancy production term within the TKE equation was based on the 
standard eddy-diffusivity formulation for the virtual potential temperature flux. This 
formulation, however, is only valid in a cloud-free environment. If a model uses this formulation 
when a dry-adiabatic layer first saturates, the virtual potential temperature gradient will be close 
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to zero and produce no additional source of TKE to increase the eddy diffusivity and mix the 
dry-adiabatic profile towards a moist-adiabatic profile. 

4.1.2.3 Improvements for Saturated Layers 

The equation for TKE in turbulent flow includes a buoyant production/destruction term 
given by 


bp - — we' 

0 V v 


(4.3) 


where g is gravity, 0 v is virtual potential temperature, and w'0' is the buoyancy flux. When 


buoyant production (BP) is positive it produces TKE and when it is negative it is a sink of TKE. 


Although it is possible to compute BP using an eddy diffusivity formulation for the flux, 


Wd v =-K 


d& v 

dz 


(4.4) 


0 v is not a conservative variable in saturated conditions. Thus, this eddy diffusivity formulation 
for the buoyancy flux is not generally valid. To obtain a general expression for buoyancy flux, 
we can use the definition of 6> v , and relate it to the conservative variables given in (4.1) and 

(4.2). Following Deardorff (1976), 


WG' V = A W9j + B Wq\ 


(4.5) 
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where all variables are as defined above, and A and B are functions of the basic thermodynamic 
variables. We can then re-write (4.5) using eddy diffusivity K as 


w'0' v = -KA ^ - KB . (4.6) 

Note that (4.6) is equivalent to (4.4) for unsaturated conditions, but for these “dry” conditions 
there are different functions for A and B. 

Because the TKE is often defined at the boundaries of the grid volumes where the 
temperature and moisture fields are defined, a problem occurs when computing the buoyancy 
production at a vertical boundary between saturated and unsaturated air. Neither the saturated or 
unsaturated formulation for buoyant production of TKE is strictly correct in this boundary region 
(e.g., Sommeria and Deardorff 1977, Munoz et al. 2000). For the purpose of this study we 
simply use the unsaturated formulation in the boundary zones. 

We also re-formulated the dissipation of TKE used by Ballard et al. (1991) to include 
moist effects in the stability parameter, N 2 , the Brunt-Vaisala frequency, in saturated layers. 
The computation of the Blackadar length scale used in the dissipation was also re-designed to 
better represent the vertical turbulent scales in saturated and unsaturated layers, based on the 
predicted TKE profile and local minima. 

We have benchmarked this TKE sub-model and saturation effects on turbulence against 
that of Duynkerke and Driedonks (1987) for two cloud-topped PBL cases. One case was forced 
by only surface buoyancy fluxes while the other case was forced by radiative forcing at the top 
of the cloud. The 1-D MM5 results compared favorably with their one-dimensional results (not 
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shown). In this section we extend our work to three dimensions to demonstrate the importance 
of saturation effects on turbulence within a mesoscale model. 

Another limitation of the original version of the TKE sub-model described by Shafran et 
al. (2000) is that cloud ice cannot be treated properly using liquid water potential temperature 
and total water defined by (4.1) and (4.2). Therefore, we have recently replaced the liquid water 
potential temperature variable in the scheme by an ice-liquid water potential temperature e a , as 

described by Tripoli and Cotton (1981), 


*n=*U 


max(7\253) 


Av j(± j 

max(7\253) 


(4.7) 


where all symbols are as above, <7, is the cloud ice mixing ratio, L, v is the latent heat for lce- 
vapor phase changes, and the max(T,253) is an empirical correction which improves diagnosis of 
6 from d d at temperatures colder than 253 K. The total water mixing ratio, modified to include 

the cloud ice mixing ratio <7, , is given by 


<7, =<?v+4f +4, • 


(4.8) 


Since an equilibrium condition cannot be assumed for the ice phase, an additional mixing 
variable was needed for the transformation of total water mixing ratio back to the native MM5 
moisture variables. Currently, the TKE submodel only allows “simple ice” (no mixed phase) 
microphysics, and a mixing equation for the cloud liquid/ice is included in the scheme, in 
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addition to the one for the total water defined by (4.8). The vapor mixing ratio is computed as a 
residual by subtracting the cloud liquid/ice mixing ratio from the total water mixing ratio. 

Additional work has focused on making the TKE sub-model more robust, as well as more 
efficient. The counter-gradient terms for the TKE flux proposed by Ballard et al. 1991 were 
found to be a source of numerical instability. These correction terms were designed to account 
for the counter-gradient TKE flux observed in the lower third of the convective PBL. However, 
after several failed attempts to make these terms numerically stable, we have decided to remove 
them from the TKE flux formulation. These counter-gradient correction terms to the TKE flux, 
which are often ignored in other TKE formulations (e.g., Musson-Genon 1987), did not appear to 
have any significant positive effect on the model solutions. 

An improved algorithm to define a smaller time step within the TKE sub-model to assure 
numerical stability has also been implemented. The surface-layer (lowest model layer) equations 
are solved explicitly while all layers above the surface layer are solved implicitly. Thus the most 
restrictive numerical stability criterion for the surface-layer diffusive equation for each north- 
south “j-slice” is used to determine the sub-model time step for that slice. This allows shallower 
surface layers to be used without needing to reduce the advection time step used to integrate the 
rest of the model equations. 

4.1.3. Results with the Improved TKE Formulations for Saturation 

We have benchmarked this new version of the TKE sub-model against that of Duynkerke 
and Driedonks (1987) for two cloud-topped PBL cases using the 1-D MM5. One case was 
forced by only surface buoyancy fluxes while the other case was forced by radiative forcing at 
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the top of the cloud. The results compared favorably with those published by Duynkerke and 
Driedonks (1987). Results from these preliminary evaluations are not shown in this report. 

An example from a 3-D MM5 simulation using the improved formulation for in-cloud 
turbulence is shown in Figure 4.1.1, which shows cloud water mixing ratio and TKE for a 
vertical cross section from San Clemente to Point Magu, off the west coast of the U.S. in August 
of 1993. The marine stratus layer shown in Figure 4. 1.1. a is associated with a region of 
enhanced TKE in Figure 4.1.1b. The larger in-cloud turbulence (compared to a similar 
experiment without the improved formulation, not shown) produces greater in-cloud mixing 
which affects the amount and distribution of the grid-resolved cloud, which in turn interacts with 
the PBL via the radiative flux divergence near the top of the cloud. 

MM5 results using the Blackadar PBL are compared with those using the new TKE sub- 
model for a sounding location midway along this cross section (Figure 4.1.2). Both model 
simulations used identical initial conditions and boundary conditions. The Blackadar PBL 
solution (Figure 4.1.2a) shows an absolutely unstable profile near the base of the cloudy region 
within the PBL, while the TKE simulation (Figure 4.1.2b) shows a more realistic moist-adiabatic 
profile through the cloud layer as expected. 

4.1.4. Discussion 

Significant errors can occur in model simulations which do not properly account for 
turbulence in saturated layers. The PSU TKE-predicting sub-model for turbulence, now 
available in MM5V3, has been reformulated and used to demonstrate the sensitivity of the model 
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1. 3-D MM5 TKE-sub-model simulation for a vertical cross section from San 
Clemente to Point Magu, California, off the west coast of the U.S., at 6 h (0600 
UTC 25 August 1993). (a) cloud-water mixing ratio (contour interval of 0.05 
g kg' 1 ) and (b) TKE (contour interval of 0.1 J kg' 1 ). 
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a. 




Figure 4.1.2. MM5 simulated sounding located midway along the cross section shown in 

Figure 4.1.1 at 6 h (0600 UTC 25 August 1993) using (a) Blackadar PBL and 
(b)PSUTKE sub-model. 
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solution to proper treatment of in-cloud turbulence, and its complex interaction with the moist 
microphysics and radiative forcing associated with a saturated layer. The recent improvements 
to the scheme, beyond that described in Shafran et al. (2000), include better representation of 
buoyant production and dissipation of TKE in saturated layers, including layers saturated with 
respect to ice. The mixing of ice-liquid potential temperature, total water mixing ratio and cloud 
water/ice mixing ratio conserves water mass and energy during vapor-liquid-ice phase changes. 

In a future release, the TKE sub-model will be adapted to mixed-phase microphysics, 
where the sum of the vapor and cloud liquid mixing ratios will be computed as a residual from 
the total water mixing ratio minus the cloud ice mixing ratio. The individual vapor and cloud 
liquid mixing ratios will then be recovered using an equilibrium assumption. 


4.2 Shallow Convection in the 3-D MM5 
4.2.1 Introductory Remarks 

Shallow clouds exert strong influences on the radiation budget of the earth- atmospheric 
system. They also play an important dynamic role in the transport of heat and moisture from the 
well-mixed planetary boundary layer to the free troposphere. In addition, by modifying the 
thermal and moisture structure of their environment, shallow clouds may affect the onset of deep 
convection and the subsequent development of mesoscale circulations. Furthermore, convective 
clouds may have strong influences on air pollution by causing boundary layer pollutants to be 
vented into higher levels, by altering chemical reaction rates through aqueous processes, and by 
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strongly affecting the actinic radiation flux. It is important that the influences of shallow clouds 
are included in the numerical models. However, most 3-D meteorological models used for air- 
quality and other mesoscale applications have treated shallow convective processes in inadequate 
ways. Often, shallow convection is represented as a simple mixing process, with cloud fraction 
diagnosed directly from the environmental relative humidity. While this approach is fairly 
reliable for global models, it is not accurate for regional models, especially for the applications 
associated with air-quality studies which require more detailed information on cloud properties 
such as cloud-base mass flux, number of clouds, cloud depth, entrainment and detrainment 
profiles. 

In order to more accurately represent the effects of shallow clouds in mesoscale 
numerical models, a shallow-convection parameterization scheme has been developed at Penn 
State University (PSU). The scheme is closely associated with boundary layer turbulent 
processes and can transition to either a deep convection scheme (Kain-Fritsch) in convectively 
unstable environments or to an explicit moisture scheme in moist stable environments. The 
scheme’s closure assumption uses a hybrid formulation based on boundary layer turbulent 
kinetic energy (TKE) and convective available potential energy (CAPE), while its convective 
trigger is primarily a function of boundary layer TKE. The scheme determines updraft 
characteristics using the vertical velocity calculated from the parcel buoyancy equation for 
entraining/detraining convective clouds. Using the updraft detrainment profile as a source term, 
neutrally buoyant cloud (NBC) is predicted by two governing equations (one for cloud fraction 
and one for liquid water content), with dissipation processes including precipitation, horizontal 
and vertical turbulent mixing, ice settling and cloud-top entrainment instability. A detailed 
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description of the shallow-convection scheme is given in Deng et al. (2002a) and in Section 2 of 
this report. 

Following earlier evaluations in a 1-D version of the Penn State/NCAR MM5 model 
(Deng et al. 2002b), the shallow-convection scheme has been installed and tested against 
observed data sets in marine environments in the 3-D MM5 (Deng et al. 2000). The evaluations 
over the mid-Atlantic convective environment using the observational data from the ASTEX has 
shown that the scheme is able to reproduce the shallow fog-stratus of mid-latitudes and the 
subtropical trade-cumulus environment over the ocean (Deng, et al. 2000). This section of the 
final report focuses on evaluating the shallow-convection scheme in a continental stratocumulus 
environment over the Southern Great Plains (SGP) region. 

4.2.2 Model Description 

Version MM5v3 of the 3-D non-hydrostatic Penn State/National Center for Atmospheric 
Research (PSU/NCAR) mesoscale model is used for these experiments. The model is configured 
with the PSU shallow-convection parameterization and the 1.5-order turbulence sub-model that 
includes saturation effects in the prediction of turbulent kinetic energy (TKE) (Stauffer et al. 
1999, Shafran et al. 2000). The turbulence scheme also is used to determine the planetary 
boundary layer (PBL) depth. Additional model physics include an explicit moisture scheme that 
predicts resolved-scale cloud and precipitation processes for water and ice (Dudhia 1989). For 
radiation, a two-stream broad-band atmospheric radiation scheme that calculates radiation fluxes 
in the column (Dudhia 1989) is used. As explained in Section 2, use of the shallow-convection 
scheme requires the column radiation subroutine to be called three times at a given time step 
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(every 15 minutes) to account for sub-grid portions of the column having clear or cloudy air. 
This is intended to calculate the radiative fluxes more accurately by considering partial 
cloudiness in the grid column. 

For these experiments, the MM5 was configured with nested domains of 108-km (37x43 
grid points) and 36-km (73x73 grid points), both centered at 35°N and 100°W. Vertically, we 
used 30 vertical layers, with 80-m layers from the surface to about 1400 m AGL. The layer 
thickness gradually increased above this height. The top of the model's reference-state pressure 
field was fixed at 100 hPa. Initial and boundary conditions were specified from NCEP global 
spectral model analysis enhanced by surface and rawinsonde data. 

4.2.3 Methodology and Case Description 

In this part of the study, we chose the case of 12-13 April 1997. Figure 4.2.1 shows 
MM5 initial sea-level pressure (hPa) at 1200 UTC, 12 April 1997 on the 36-km domain. Figure 
4.2.2 shows a visible satellite image of stratocumulus clouds over the Southern Great Plains 
(SGP) region at 1709UTC, 12 April, centered on the ARM-CART central facility (36.61°N, 
97.49°W). The entire domain of the satellite image covers the rectangular region in Fig. 4.2.1. 
A strong cold front had passed the central facility around 0000UTC, 1 1 April, bringing north- 
northwesterly winds and sharply falling temperatures. By early morning, at 1200 UTC, the 
surface temperature dropped from -20 °C on the previous afternoon to -3 °C at the central 
facility, and post-frontal stratus clouds covered the region. Cold advection and north- 
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Figure 4.2.1. Analysis of sea-level pressure (hPa) at 0000 UTC, 12 April 1997 on the 36-km 
domain. The rectangular area over Kansas-Oklahoma region corresponds to the area covered by 
the observed satellite image in Fig. 4.2.2. Dashed line in the rectangular area shows the location 
for the cross sections in Fig. 4.2.6. The open circle on the dashed line indicates the location 
where vertical profiles in Fig. 4.2.7 are created. 




Figure 4.2.2. Visible satellite image at 1709 UTC, 12 April 1997, showing widespread 
stratocumulus clouds over the Kansas-Oklahoma region, and decreasing cloud cover towards the 
southwest. 
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northwesterly winds prevailed over the site through the next 36 h as high surface pressures were 
established over KS and OK. Following sunrise, heating of the surface led to the growth of the 
PBL and lifting of the cloud base. Extensive stratocumulus clouds evolved during the period in 
response to the diurnal heating cycle and boundary layer forcing. Figure 4.2.3 shows the 
observed cloud base, cloud top and PBL depth at the central facility, measured using a laser 
ceilometer, micropulse lidar and 35-GHz cloud radar. The cloud layer was initially very low, 
then gradually lifted as the boundary layer grew. As dry air was entrained from above into the 
growing PBL, the clouds changed from initial low stratus to broken stratocumulus clouds that 
extended across most of the region. These stratus and stratocumulus clouds accompanying the 
cold-air outbreak (Fig. 4.2.2) persisted through the day (12 April) and then dissipated near sunset 
as the surface heat flux became negative and a stable boundary layer replaced the daytime 
turbulent mixed layer. 



Figure 4.2.3. Observed cloud base (dots, m), cloud top (crosses, m) and mixed-layer depth 
(open circles, m) from 1200UTC, 12 April to 1200 UTC, 13 April 1997, at the Southern Great 
Plains ARM-CART central facility (36.61°N, 97.49°W). 
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To demonstrate the performance of the PSU shallow-convection scheme, we focus on the 
36-km domain and conduct two 3-D experiments on this domain. The first (Exp. 4.2.1) 
experiment is run without the shallow-convection sub-model, so that all clouds must be 
represented on the resolved scale using MM5’s explicit microphysics scheme (in this case a 
simple water/ice parameterization is used). The second (Exp. 4.2.2) experiment adds the 
shallow-convection sub-model, so that the sub-grid scale distribution of fractional cloud area and 
water/ice content are predicted. Both experiments are initialized at 0000UTC, 12 April and run 
for 24 hours. Simulated results are discussed in the following section. 

4.2.4 Results for the Stratocumulus Case 

The following discussions will focus on results after sunrise between 1200 UTC, 12 April 
and 0000 UTC, 13 April on the 36-km domain. Figure 4.2.4a and 4.2.4b show the resolved scale 
liquid water content (scaled by 10000) of low clouds in Exp. 4.2.1 at 1800 UTC (local noon) and 
2100 UTC (mid-afternoon), respectively. The first contour, 0.005 g kg' 1 gives the approximate 
outline of the low stratus for Exp. 4.2.1. Since there is no mechanism in Exp. 4.2.1 for 
producing partial cloudiness, Fig. 4.2.4a and 4.2.4b show solid stratus clouds with low water 
content across northern KS, but virtually no clouds over OK. That is, there are only two options: 
completely cloudy or completely clear. While the model captures the general observed trend of 
less cloudiness to the south and more to the north over the SGP, the true character of this case 
cannot be represented in the model (without a shallow-convection parameterization). 
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Figure 4.2.4. MM5-predicted liquid/ice water content (g kg' 1 ) of explicit cloud (scaled by 
10000) at 1 km AGL in Exp. 4.2.1: a). 1800 UTC. b). 2100 UTC, 12 April 1997. Contour 
interval is 0.005 g kg' 1 . Rectangular area corresponds to the area covered by the satellite image 
in Fig. 4.2.2. 
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Next, Fig. 4.2.5a and 4.2.5b give the fractional cloud area (%) for low clouds predicted in 
Exp. 4.2.2, with the PSU shallow-convection scheme, at 1800 UTC and 2100 UTC, 12 April 
1997. Contour interval is 10%. Comparing Figs. 4.2.4 and 4.2.5, we Find that the difference is 
dramatic. Both experiments show 100% cloud cover in northeastern KS gradually trending 
toward fewer clouds to the south, with clear skies over central TX (beginning at the southwest 
comer of the box indicating the domain of the satellite image in Fig. 4.2.2). Comparison with 
Fig. 4.2.2 confirms that this pattern is correct. Although Fig. 4.2.5b is valid about 4 hours after 
the satellite picture, it shows that the model generates the strongest gradient in the cloud fraction 
over the central part of the verification region. 
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Figure 4.2.5. MM5-predicted fractional cloud area (%) for low clouds in Exp. 4.2.2: a). 1800 
UTC. b). 2100 UTC, 12 April 1997. Contour interval is 10%. Rectangular area corresponds to 
the area covered by the satellite image in Fig. 4.2.2. 
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In order to examine the evolution of the 3-D cloud fields, we make a northwest-southeast 


cross section along the dashed line in the box covered by the satellite picture. Figure 4.2.6a, b, c 
and d show the cross sections of fractional cloud area (%) in Exp. 4.2.2 at I500UTC, 1800UTC, 
2100UTC, 12 April and 0000UTC, 13 April, respectively. Initially at 1500UTC (shortly after 
sunrise. Fig. 4.2.6a), solid stratus dominates the entire region. By 1800UTC (local noon) the 
stratus layer has broken and the cloud layer has risen significantly. This is consistent with the 
observation shown in Fig. 4.2.3. The maximum cloud fraction is less than 50%. By mid- 
afternoon (Fig. 4.2.6c), as updrafts continue to detrain their mass into their environment to form 
NBC, the maximum cloud fraction gradually increased to 85-90% that is consistent with the 
observed stratocumulus clouds in the satellite images shown in Fig. 4.2.2. Around sunset at 
0000 UTC, 13 April, stratocumulus clouds start to dissipate, and eventually disappear over the 
entire region as the PBL collapses. 


162 





Figure 4.2.6. MM5 north-south cross section of predicted fractional cloud area (%) versus 
pressure (hPa) in Exp. 4.2.2 along the dashed line in Figure 4.2.1: a). 1500UTC, b) 1800UTC, c) 
2100 UTC, 12 April 1997, d) 0000 UTC, 13 April 1997. Contour interval is 10%. Dashed curve 
is the MM5 predicted PBL depth. 
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The behavior of the shallow-convection scheme is further examined by looking at vertical 
profiles of cloud fraction and cloud water calculated using a prognostic scheme (Deng et al. 
2002a). Figure 4.2.7a shows the vertical profile of cloud fraction at location on the cross section 
described above (shown as open circle on the dashed line in Fig. 4.2.1) at 1800 UTC. At this 
time the maximum NBC fraction is about 50% although the environment relative humidity is 
about 90%, with cloud top at 1570m and LCL at 634m. It is also shown in this figure that the 
fraction of updraft is only less than 5%. In the vertical profile of cloud water/ice content (Fig. 
4.2.7b) at this time, the maximum cloud water content in the updraft is about 0.9 g kg' 1 (shown 
as L u ), while the water/ice content in the NBC is less than 0.01 g kg' 1 shown as L c ). These 
values are expected in this type of cold-season stratocumulus clouds. 

4.2.5 Discussion 

The PSU shallow-convection scheme has been tested in a stratocumulus environment 
over the Southern Great Plains (SGP) region during a case of 12-13 April 1997. Comparison 
between simulations with and without the shallow-convection scheme was performed. With the 
shallow-convection scheme, the MM5 model was able to reproduce the observed characteristics 
of the stratocumulus environment over land. The simulated results including the cloud fraction, 
cloud water/ice content were generally consistent with the observations from the field 
experiments (including radar and satellite imageries). 
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Figure 4.2.7 MM5 predicted cloud vertical profiles in Exp. 4.2.2 at 1800 UTC, 12 April 1997 
at location marked by an open circle on the dashed line in Fig. 4.2.1. a). Cloud fraction (%) and 
resolved-scale relative humidity (%), b). Cloud water/ice content (g kg' 1 ). 




4.3 


PLACE and Shallow Convection in the 3-D MM5 


4.3.1 Overview of the System 

The basis of the MM5-PLACE development at Penn State is version 3 of the MM5 
modeling system (Dudhia 1993, Grell et al. 1994). This version has the capability of using the 
Oregon State University (OSU) land-surface model, so that its preprocessor programs can 
produce several types of data that are required by PLACE or other complex land-surface models. 
The implementation of PLACE into MM5 has been designed so as to minimize the number of 
changes necessary in the standard public-release version of the MM5 modeling system. The 
input data that are specific to PLACE are ingested through an input file with the same data 
format as the standard MM5 input files, but which has appended all required PLACE variables. 
For example, these variables include the initial soil moisture and temperature fields that come 
from the offline version of PLACE described in Section 3.2. 

When PLACE is not used, the simpler Blackadar force-restore Slab model (Zhang and 
Anthes 1982) is used to predict ground temperature. In the Slab land-surface model, the ground 
temperature of a thin “slab” representing both soil and vegetation effects changes in time 
according to a surface energy balance. This predictive ground-temperature equation accounts for 
differences in surface characteristics (e.g., roughness, albedo, moisture availability) as 
represented through a lookup table based on climatology, and the cloud and radiation effects 
computed by the MM5 explicit moisture and column radiation submodels. The surface moisture 
availability is specified at each grid cell and it does not change with time as the soil moisture 
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content profile does in PLACE. Thus, the surface moisture in slab does not respond to model- 
predicted precipitation or the lack thereof as it does in PLACE. 

The 3-D MM5 system is used in this section to investigate the interactions between 
surface processes, the boundary layer and shallow clouds. In addition to the PLACE land- 
surface model, study of these interactions depends on two additional physical parameterizations: 
the Penn State 1. 5-Order TKE turbulence scheme (Stauffer et al. 1999, Shafran et al. 2000), and 
the Penn State shallow-convection parameterization developed by Deng et al. (2002a). The TKE 
scheme and the improvements introduced during this study were described in Section 4.1. The 
shallow-convection scheme was described fully in Section 2, and improvements to the original 
scheme developed under this project were described in detail in Section 3.3. 

The shallow-convective scheme can be considered an extension of the Kain-Fritsch (K-F) 
deep-convective parameterization, for cases of non-precipitating clouds. The shallow- 
convection scheme considers the effects of sub-grid clouds that are too shallow to activate the 
standard K-F scheme. Shallow clouds are triggered when parcels having thermodynamic 
properties defined in the near-surface layers are able to reach their lifting condensation level. 

The closure algorithm that determines the magnitude of the cloud-base mass flux distinguishes 
between a shallow, a hybrid, and a deep convection regime. In the purely shallow regime the 
mass flux is related to the TKE in the boundary layer (BL). In the deep convection regime the 
CAPE-removal closure of K-F is adopted, while in the hybrid regime an interpolation between 
the two extreme regimes is used. The cloud model used to calculate the cloud-environment 
interaction is similar to the one used in K-F. Finally, the shallow-convection scheme allows for 
sub-grid cloud water to be detrained from the active updrafts computed by the cloud model into 
sub-grid neutrally buoyant clouds (NBCs). Prognostic equations for the evolution of the sub-grid 
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NBC water content and cloud fraction are solved. These equations include the effects of 
evaporation, turbulent mixing, drizzle, etc. These sub-grid cloud fields in turn are also 
considered by the radiative parameterization and therefore can affect the energy input at the 
surface. 


4.3.2 Case Description 

In this section we show the results from the coupled MM5v3-PLACE with shallow 
convection for the day of 6 July 1997, focusing on the Central Great Plains region. Figure 4.3.1 
shows the sea-level pressure analysis at the initial time (12 UTC) on the 36-km MM5 domain. 
Synoptic pressure gradients over the Central Plains were weak on this day, with a trough 
extending from the Great Lakes to northeastern Kansas. Through the night deep convection 
developed over northeastern Texas and spread into southeastern Oklahoma during the day. The 
rest of Oklahoma developed a field of shallow convection on 6 July with no precipitation. This 
spatial gradient of convective conditions represents an interesting and challenging environment 
for testing and improving the parameterizations used in this work. 

4.3.3 Experimental Design 

The numerical experiments described here are designed to assess how these different 
physical processes interact and affect the dynamics of the BL. There are four basic runs for 
different combinations of the surface and convective parameterizations (runs R04, R06, R13, and 
R14). Run R04 uses PLACE and the Kain-Fritsch (K-F) deep convection parameterization, R14 
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Figure 4.3. 1 Model domain map and sea level pressure analysis at 12 UTC 6 July 1997 
(contours in hPa). 
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uses PLACE and K-F with the Penn State shallow-convective parameterization. Run R06 uses 
the slab force-restore model and the K-F deep convective scheme, and R13 uses slab with the K- 
F deep convective scheme and the shallow-convective scheme. There are two additional 
sensitivity runs. Run R15 is similar to R04, but the model precipitation is not allowed to affect 
the soil moisture. Run R17 is similar to R13, but the shallow cloud-base mass flux has been 
augmented by a factor of three. All runs have been made with a 36-km grid and 32 levels in the 
vertical extending up to 50 hPa and with 15 levels below 850 hPa. All runs are integrated from 
12 UTC 6 July 1997 to 0000 UTC 7 July 1997. 

4.3.4 Results 

4.3.4.1 Shallow Cloud Fields and Their Effects 

Figure 4.3.2 shows the model's low-level sub-grid cloud area fractions predicted at 18 
UTC July 6, 1997, 6 hours after initialization. The plotted field corresponds to the maximum 
cloud fraction predicted by the model for the pressure levels between 970 and 800 hPa. Figure 
4.3.3 shows the low cloud amounts derived from GOES satellite data by Patrick Minnis at 
NASA Langley Research Center (data obtained from http: //www-pm. larc .nasa . gov/). 
Comparison of Figs. 4.3.2 and 4.3.3 shows that the spatial distribution of low-level cloudiness 
has been partially captured by the parameterization. The SW-NE band that crosses Oklahoma is 
reproduced in the model, although shifted eastwards somewhat. The predicted cloud amounts, 
however, are much smaller than those derived from satellite. The sensitivity run R17 (with an 
augmented cloud mass flux) shows the same spatial pattern, but peak amounts around 40% 
compared to the peaks around 30% shown in Fig. 4.3.2. The runs that used PLACE as the land- 
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Run 1 3 Shallow Low Cloud Fractlon(%) 970706 18UTC 
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Figure 4.3.2 Model calculated sub-grid low-level cloud fractions in run R13 for 6 July 1997 at 
18 UTC. Run R13 uses SLAB surface model and Kain-Fritsch convective parameterization. 
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Figure 4.3.3 
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GOES derived low level cloud fractions for 6 July 1997 at approx. 18:24 UTC. 
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RUN06 970706 1 8UTC SWD vs MaxCLW RUN06 970706 1 8UTC SWD vs MaxCLW 
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Figure 4.3.4 Scatter plots of model downward shortwave radiation fluxes at the surface for 6 
July 1997 at 18 UTC. a) Abscissa variable is column maximum resolved-scale cloud water 
content in run using KF convective scheme only, b) as a) but for grid points with maximum 
cloud water content below 10' 5 kg/kg, c) as b) but for run that uses shallow convective scheme, 
d) as c) but abscissa variable is sub-grid cloud area fraction. 
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surface parameterization showed a more extensive region of shallow clouds, associated with 
larger evaporation rates computed in those runs (not shown). 

One of the physical processes that a sub-grid cloudiness scheme tries to represent is to 
account for the effect of unresolved clouds on the solar radiation reaching the surface. Figure 
4.3.4 shows the dependence of model downward solar radiation (SWDOWN) as a function of 
resolved and sub-grid cloud variables, for all grid points at 18 UTC. Panel a shows SWDOWN 
versus the maximum resolved cloud water content in the column, for a run without the shallow 
convection scheme. As expected, larger cloud water amounts are associated with reduced 
radiative fluxes at the surface. Panel b is a blowup of the same plot, but for maximum cloud 
water amounts lower than 10' 5 kg/kg. In these clear and almost-clear conditions the standard 
model computes radiative fluxes that are consistently large, since there is no sub-grid cloudiness 
to produce attenuation. Panel c is constructed as panel b, but for results of run R13, that includes 
the use of the shallow-convection parameterization. In this case there is more variability in the 
radiative fluxes for the grid points that are almost clear with respect to resolved cloud water. 
Finally, panel d uses the same quasi-clear points of panel c, but the scatter is between SWDOWN 
and sub-grid cloud area fraction. The latter variable explains much of the computed attenuation 
of the solar radiation reaching the surface. 

The effects of the shallow-convection parameterization on the vertical redistribution of 
water vapor and temperature are illustrated by the vertical profiles shown in Fig. 4.3.5. Panels a) 
and b) show the thermodynamic profiles for point B5 in east-central Oklahoma, and panels c) 
and d) are for point B6 in central-southern Oklahoma (locations shown in Fig. 4.3.2). Observed 
soundings are plotted in each panel as bold lines. Especially over point B5, where shallow 
clouds were more prevalent in the model, the shallow convection parameterization appears to dry 
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Figure 4.3.5 Observed and modeled vertical profiles of thermodynamics variables, a) potential 
temperature at point B5, b) water vapor mixing ratio at point B5, c) as a) but for point B6, d) as 
b) but for point B6. Results are for runs with different convective schemes. 
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the upper portion of the boundary layer in R13 and moisten the layers above it. With the cloud- 
base mass flux artificially increased in R17, the drying effect extends throughout the BL. The 
temperature profile is affected by the vertical transport and also the radiative effects described 
earlier, and the net result compared to B5 appears to be a cooling of the boundary layer and the 
air above it. In these three-dimensional runs advection can also be important, so that these 
preliminary effects have to be studied further. 


4.3.4.2 Effects of Surface Processes 

In this section we describe the effects of using different treatments of the surface 
boundary conditions. As discussed in Section 4.3.1, the standard MM5 model uses the Slab 
parameterization to calculate the surface fluxes, including a prognostic equation for the soil 
surface temperature (only one layer) and a fixed climatological moisture availability field to 
calculate the latent heat fluxes. The alternative scheme uses the PLACE land-surface model, 
which includes tendency equations for five soil moisture layers, and seven soil temperature 
layers. 

During this time of the year Oklahoma has a strong longitudinal gradient in vegetation 
cover and soil moisture. Figure 4.3.6b shows longitudinal averages of water vapor mixing ratios 
derived from surface Oklahoma Mesonet stations. Each point is constructed by averaging the 
measurements of all stations located in a 1 degree latitude by 1 degree longitude box centered at 
36° N and moving -99.5° W to -95.0° W. The measured vapor gradients at 12 UTC and 18 UTC 
show that the eastern part of the state develops a moister boundary layer at 1 8 UTC than the 
central and western parts, consistent with the vegetation cover gradient. Model results shown in 
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Figure 4.3.6 Longitudinal cross sections of surface layer water vapor at 36 N for 6 July 1997 at 
12 UTC and 18 UTC. a) Model results at first layer, b) Averages derived from MESONET 
observations. 
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Fig. 4.3.6a indicate that the standard slab surface model does not capture this gradient in 
moisture (run R06). The run R04 that uses PLACE does indeed produce more moistening in the 
eastern part of the state, although the mixing ratios appear to be too high. In Run R04, the K-F 
deep-convection scheme produced widespread light rainfall over Oklahoma (largely due to 
inaccurate initialization of the precipitation already occurring at the intial time, 12 UTC), which 
caused an increase of the surface soil moisture and subsequent high evaporation rates calculated 
by PLACE. This contributed to the excessive boundary-layer water vapor in R04 at 18 UTC. 

Run R15 does not include the rain feedback on the soil moisture (and therefore on the 
evaporation rates). Thus, experiment R 15 also shows the longitudinal gradient of moisture, but 
with lower overall values of mixing ratio and perhaps too much drying in the central region of 
the state. Dynamic initialization can be used in the future to minimize these initialization errors. 

The diurnal cycle of water vapor mixing ratio in the convective boundary layer usually 
shows a moistening in the early morning hours and a subsequent drying due to entrainment of the 
drier air from above the BL. If the moisture supply at the surface is high, this later drying can be 
reduced. Figure 4.3.7 shows in bold lines the surface temperature and mixing ratio measured at 
the Mesonet stations closest to points B5 (panels a and b) and B6 (panels c and d). The second 
station near B6 shows a more pronounced drying in the afternoon than the first station near B5 
does. From the model runs, the run R15 that uses the PLACE model with no rain feedback 
appears to produce the best time series of temperature and mixing ratio over both points. 

Figure 4.3.8 shows potential temperature and water vapor mixing ratio vertical profiles at 
18 UTC 6 July 1997, over points B5 (panels a and b) and B6 (panels c and d), measured and 
calculated by the model runs that used the PLACE land-surface sub-model. 
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Figure 4.3.7 Time series of temperature and water vapor mixing ratio, a) Mesonet station and 
model grid point closest to point B5, b) as a) but for water vapor mixing ratio, c) Mesonet station 
and model grid point closes to point B6, d) as c) but for water vapor mixing ratio. 
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Figure 4.3.8 Observed and modeled vertical profiles of thermodynamics variables, a) potential 
temperature at point B5, b) water vapor mixing ratio at point B5, c) as a) but for point B6, d) as 
b) but for point B6. Results for runs with different surface schemes. 
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Comparison of these soundings with those in Fig. 4.3.5 shows a tendency for the boundary layer 
calculated by PLACE to be moister and cooler than the observations. The bias is reduced when 
the rain feedback on the evaporation is eliminated (run R15). The measured temperature profile 
over point B6 has the peculiarity of having a strong unstable surface layer, but a stably stratified 
and not well-mixed boundary layer. It is possible that this condition is due to advection that is 
not captured adequately by the model results. 

The best way to test the performance of the surface parameterizations would be by 
comparing predicted surface fluxes with observations. However, reliable measurements of 
surface fluxes are scarce. For the period of interest the SGP data archive provides only 8 flux 
measurement points with no quality flags. Figure 4.3.9 shows modeled and measured surface- 
flux time series for points E09 (panels a and b) and E16 (panels c and d) (locations shown in Fig. 
4.3.2). In the region where these and the rest of the measurements are available, the model using 
PLACE tends to underestimate the evaporation rates, as in the case of the results for station E09. 
The comparison of the measurements with the PLACE runs is better at station E16. However at 
these points we see a tendency for the modeled latent fluxes to peak earlier in the morning and 
then decrease, a feature that has to be better understood with more detailed analysis of the 
results. 


4.3.4.3 Effects of Sub-grid Soil Moisture Heterogeneity 

The latent heat fluxes calculated by the 3-D MM5-PLACE model in the previous section 
show some features that we will try to describe and explain here in more detail. The time series 
of modeled latent heat fluxes in panels a) and c) of Fig. 4.3.9 show an early peak in evaporation 
and then a tendency for near-constant or declining values in the afternoon. Although in some 
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Model Heat Fluxes at E09 EBBR Measured Heat Fluxes at E09 



Model Heat Fluxes at El 6 ECOR Measured Heat Fluxes at E16 



Figure 4.3.9 Model and observed surface heat fluxes, a) Model sensible and latent heat fluxes 
at station E09, b) observed fluxes at E09 with Energy Balance-Bowen Ratio (EBBR) instrument, 
c) model sensible and latent heat fluxes at station E16, d) observed fluxes at E16 with eddy 
correlation instrument. See Fig. 4.3.2 for location of stations. 
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cases the observations show a similar behavior (panel d of the same figure), in many other cases 
the diurnal cycle of evaporation is more symmetric with respect to noon (panel b). Also, in 
visualizing horizontal maps of calculated evaporation rates, there is a tendency to see sharp 
gradients between zones of little evaporation and other zones in which evaporation rates are 
large. We try to understand these model results and explore the effect that sub-grid soil moisture 
variability can have in this respect. 

To illustrate two different types of model results, we present in Fig. 4.3.10 modeled latent 
heat fluxes for two grid points, those closest to ARM-CART stations Cl and B5. Panels a) and 
c) present results for Run R15 that did not include precipitation feedback into the soil moisture, 
and also did not include sub-grid soil moisture variability. In each panel we present the total 
evaporation rate, as well as the parts due to evaporation from the soil (very small in both cases) 
and vegetation. Panels b) and d) show the same results, but for Run 04s4 that includes 
consideration of sub-grid soil moisture. The discussion of how we estimated soil moisture 
variability is presented in the following section. For now we concentrate on panels a) and c) and 
try to explain the significant difference in the diumal cycle of evaporation in these two grid 
points. In particular we would like to understand for Cl (panel a) its early peak in evaporation 
rate and the leveling off in the afternoon, and for B5 (panel c) its well defined diumal cycle and 
higher values of evaporation rates. 

The different time traces of model evaporation at points Cl and B5 can be due to 
differences in the soil and vegetation types (and associated parameters), or different states of the 
soil variables (moisture and temperature), or a combination of both. Figures 4.3.11 and 4.3.12 
show the time evolution of soil moisture and soil temperature (upper layer) for the same model 
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Runl 5 Evaporation Point Cl Run04s4 Evaporation Point Cl 
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Figure 4.3.10 Time series of evaporation rates for 6 July 1997. a) Total latent heat flux (LHF), 
soil evaporation (Soil), and plant transpiration (Trans) at grid point closes to Cl for Run 15 that 
has no rain feedback on soil moisture and no soil moisture heterogeneity, b) as a) but for Run 
04s4 that considers soil moisture heterogeneity, c) as a) but for gridpoint closest to point B5. d) 
as b) but for gridpoint closest to point B5. 
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Figure 4.3. 1 1 Soil moisture at five soil layers for MM5-PLACE computations, a) For point 
closes to Cl and run R15 that does not include variability in soil moisture, b) as a) but for run 
R04 that includes soil moisture variability, c) as a) but for point B5, d) as b) but for point B5. 
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Figure 4.3.12 Soil temperature for MM5-PLACE computations, a) For point closest to Cl and 
run that does not include variability in soil moisture, b) as a) but for run that includes soil 
moisture variability, c) as a) but for point B5, d) as b) but for point B5. 
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runs as in Fig. 4.3.10. It is apparent that the upper soil layers are initially moister at point B5 
than at point Cl. Soil temperature reaches higher values at Cl than at B5, consistent with the 
smaller evaporation fluxes at Cl. In terms of differences of the parameters characterizing soil 
and vegetation at both grid points, both points have the same soil type (silt loam) and this 
eliminates a potentially important source of differences between them. However, point Cl has 
mixed dryland/irrigated cropland and pasture and point B5 has savanna, and the vegetation 
fraction at Cl is 0.41 and that at B5 is 0.79. 

In order to study in detail the behavior of the evaporation rates calculated by PLACE, we 
isolated its evaporation subroutine and constructed a Matlab graphical user interface (GUI) that 
allows for easy interactive computation of evaporation rates for different conditions in the 
surface parameters and forcings. In particular we calculated evaporation rates as function of soil 
temperature, soil moisture, or soil moisture variability. Figure 4.3.13 shows evaporation rates as 
a function of soil temperature for points Cl (panel a) and B5 (panel b). The bold lines indicate 
total evaporation. The dash-dotted line indicates the fraction that is due to evaporation from 
plant storage. The GUI results presented in Fig. 4.3.13a show for point Cl a pattern similar to 
that of the 3-D MM5-PLACE computations (Fig. 4.3.10a). Total evaporation peaks early and 
then flattens out to a nearly constant value. (Note that the temperature abscissa in the figure can 
be considered also a time variable, since we started the computation with the temperature and 
plant storage conditions similar to the 12 UTC conditions, shortly after the daily minimum 
temperature, and we finished with the maximum temperature observed around 18 UTC.) A 
detailed analysis of the GUI results shows that the early peak in evaporation rate is due to 
depletion of plant water storage in the model. Evaporation of plant storage in the model occurs 
when a layer of soil reaches its maximum evaporation rate. Plants are assumed to continue 
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Figure 4.3. 13 Computed evaporation rates as a function of soil temperature, a) for conditions 
similar to point Cl, b) for conditions similar to point B5. 
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evaporating at the potential rate, so that any stored water is depleted. After some time the plant 
storage vanishes and the evaporation rate returns to the maximum sustainable rate that can be 
taken by the plants directly from the soil (threshold evapotranspiration). Results for point B5 in 
panel b) show that the threshold plant evaporation is larger than for Cl. This is due to a moister 
soil and to a larger vegetation cover fraction. Consistently, the depletion of plant storage occurs 
later than for point Cl. The results in the 3-D MM5-PLACE simulation for B5 (Fig. 4.3.10c) do 
not show the peak and flat regimes of the evaporation rates because the soil temperature (Fig. 
4.3.12c) reaches a maximum of only about 301 K, and Fig. 4.3.13b shows that the peak and flat 
regime occurs at temperatures several degrees warmer. 

In summary, the diurnal cycle of evaporation in the PLACE model is seen to occur in 
three stages. In the early morning the evaporation rates increase at the potential rate (stage 1). 

At some point some of the soil levels reach their threshold evaporation rates. At this point plant 
storage begins to be depleted, so that the net evaporation rate continues to increase at the 
potential rate (stage 2). Finally, when plant storage is depleted, the total evaporation rate returns 
back to the threshold value (stage 3). 

Next, we turn to the question of how the inclusion of sub-grid soil moisture variability 
may alter the diurnal cycle of evaporation in the model. Figure 4.3.14 shows the results obtained 
with the GUI for evaporation rates as a function of soil moisture (panels a and b), and as function 
of the standard deviation of soil moisture (panels c and d). The latter are calculated for 
conditions similar as those at 18 UTC in the 3-D MM5-PLACE model results. In each panel 
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Evaporation vs Soil Moisture Cl Evaporation vs Soil Moisture B5 



Evaporation vs Std.Dev. Soil Moisture Cl Evaporation vs Std.Dev. Soil Moisture B5 



Figure 4.3.14 Computed evaporation rates as a function of soil moisture and soil moisture 
variability, a) Conditions like in point Cl, abscissa is soil moisture, b) as a) but for B5, c) 
conditions like in point B5, abscissa is standard deviation of soil moisture, d) as c) but for B5. 
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we plot the total evaporation rate as a bold line along with the fractions from soil versus 
vegetation. Soil and vegetation evaporations have a very different sensitivity to soil moisture. 
Soil evaporation increases rapidly, while vegetation evaporation increases a little bit slower. 
Depending on the shape of the curves in panels a) and b), the inclusion of soil moisture 
variability can produce an increase, a decrease or not greatly affect the evaporation rates for soil 
and vegetation. Analysis of panels c) and d) shows that for conditions similar to point Cl the 
inclusion of soil moisture variability greatly increases evaporation from soil, but it hardly 
changes the vegetation evaporation. For point B5 the effects in soil and vegetation evaporations 
have different signs and therefore the total evaporation is little affected by inclusion of soil 
moisture variability. These results explain the differences in the latent heat fluxes shown in Fig. 
4.3.10. 


4.3.4.4 Methodology to Define Subgrid Soil Moisture Heterogeneity 

Subgrid heterogeneity is addressed in PLACE either through mosaic tiling, where each 
model grid cell is divided into smaller tiles for flux calculations, and then a composite flux is 
computed for the grid cell, or a statistical specification of the variability is performed. Either 
methodology or both can be used. We use the latter statistical methodology. In order to define 
realistic subgrid soil moisture heterogeneity statistics within the 3-D MM5-PLACE model run 
that included its effects (Run04s4), we used high-resolution soil moisture data available from 
ESTAR measurements for some days of July 1997 during SGP97. Figure 4.3.15 shows scatter 
plots of the standard deviation of soil moisture versus mean soil moisture for 36 X 36 km and 12 
X 12 km averages. Based on the 36-km panel (Fig. 4.3.15a) we assigned for our 36-km MM5- 
PLACE work a maximum standard deviation of 0.10 and a coefficient of variation of 0.33 
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Figure 4.3.15 Soil moisture variability as derived from ESTAR data. Abscissa is mean soil 
moisture and ordinate is standard deviation of soil moisture, a) for 36x36 km 2 averages, b) for 
12x12 km 2 averages. 
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defined from the slope of the pattern; that is, the standard deviation divided by the mean, as 
defined in Wetzel and Chang (1988). 

The effects of including subgrid variability in the soil moisture are assessed in Figs. 
4.3.16 and 4.3.17. In the first we show histograms for latent heat fluxes at 18 UTC calculated 
over the grid points covering Kansas and Oklahoma. Panels a and b show results for total 
evaporation rates for the run without soil moisture variability and the one with it, respectively. 
When subgrid variability is not considered the histograms have a double peak structure, with 
maxima at lower and higher values of fluxes. 

Inclusion of variability reduces the peak at low values. Panels c and d shows the 
respective histograms for soil evaporation and panels e and f do the same for evaporation from 
vegetation. In both cases the frequency of very small fluxes is decreased, especially for soil 
evaporation. The corresponding histograms for sensible heat fluxes (Fig. 4.3.17) show a similar 
behavior. A double peak structure is present when no soil moisture sub-grid variability is 
considered and the histogram is more unimodal histogram as sub-grid variability is included. 
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Figure 4.3.16 Histograms of modeled latent heat fluxes, a) Results for Run 15 that does not 
include soil moisture variability, b) as a) for Runl4s4 that includes soil moisture variability, c) as 
a) but for evaporation of soil, d) as b) but for evaporation of soil, e) as c) but for evaporation 
from plants, f) as d) but for evaporation from plants. 
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Figure 4.3.17 Histograms of modeled sensible heat fluxes, a) results for Run 15 that does not 
include soil moisture variability, b) as a) for Runl4s4 that includes soil-moisture variability 
maxima at lower and higher values of fluxes. 
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5. 


SUMMARY 


The goal of this project was to develop and test a coupled model-physics suite to 
investigate the interactions between the surface fluxes, the diumal cycle, the planetary boundary 
layer characteristics and the moisture cycle (including shallow and deep clouds) in a 
nonhydrostatic three-dimensional (3-D) mesoscale model (the Penn State University / National 
Center for Atmospheric Research (PSU/NCAR) mesoscale model, widely know as MM5). 

Atmospheric modeling over continents is incomplete without adequate treatment of the 
land-surface boundary, the planetary boundary layer and its associated cloud systems. It was 
hypothesized that an improved treatment of the regional water cycle can be achieved by using a 
3-D mesoscale numerical model having high-quality parameterizations for the key physical 
processes controlling the water cycle. These included a detailed land-surface parameterization 
(the Parameterization for Land- Atmosphere-Cloud Exchange [PLACE] sub-model of Wetzel and 
Boone, 1995), an advanced boundary- layer parameterization (the 1. 5-order turbulent kinetic 
energy [TKE] predictive scheme of Shafran et al„ 2000), and a more complete shallow 
convection parameterization (the hybrid-closure scheme of Deng et al., 2002a,b) than are 
available in most current models. PLACE is a product of researchers working at NASA s 
Goddard Center in Greenbelt, MD. The TKE and shallow-convection schemes are the result of 
model development at Penn State. 

The work conducted in this study provided the first insertion of the PLACE scheme into 
the MM5v3 model (the latest version of the MM5). It was found through extensive testing in 3- 
D and 1-D frameworks that the PLACE land-surface scheme performed at about the same level 
of accuracy as other land-surface schemes used to represent the surface fluxes over the Great 
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Plains in normal operations. In some cases, PLACE was able to perform better than alternative 
schemes. However, in many situations, it was found that a more simple land-surface contained 
in the Blackadar boundary layer parameterization (Zhang and Anthes 1982), known as the 
Blackadar Slab scheme, generally provided about the same level of skill or better. Thus, while 
the results of MM5-PLACE simulations were as good or better than any others known to be done 
with PLACE, it proved difficult to show an overall clear advantage when PLACE is used. 

The underlying reason for this result is attributed mostly to the uncertainty involved in 
defining many of the key soil and vegetation parameters that must be defined in order to run 
PLACE. Even with the best soil-characteristic data sets currently available, such as the 
STATSGO data base, PLACE produced soil-moisture evolutions that did not always agree well 
with soil-moisture measurements over the Great Plains. Work is currently underway at Penn 
State to evaluate MM5-PLACE during July 1997 using ESTAR soil moisture and NDVI-based 
leaf-area index and vegetation fraction. 

The 1.5-order TKE-predicting turbulence scheme of Shafran et al. (2000) worked 
reasonably well in dry (unsaturated) conditions, but was found to have deficiencies in saturated 
layers. A reformulation was introduced that now accounts for the effects of saturation on the 
calculations in the buoyancy-production term of the turbulent kinetic energy equation, and the 
stability parameters and turbulent length scales (Stauffer et al. 1999). These improvements to 
the TKE PBL scheme were given to scientists at NASA GSFC for the collaborative work with 
Penn State reported in Lynn et al. (2001). 

This saturation dependency is missing in many current numerical models. The change 
allows the turbulence sub-model to generate more intense mixing in clouds, especially near cloud 
tops where the longwave radiative flux divergence leads can cause intense cooling. The 
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consequence of this improvement is a more realistic thermodynamic and moisture structure in 

cloudy regions and better representation of fog layers. 

Extensive evaluations and testing of the shallow-cloud parameterization also were 
conducted in the course of this investigation. It was found that this convection sub-model 
contained several layer-dependent formulations for some key quantities related to the cloud-base 
mass flux. In particular, the discrete thickness of the model layers tended to cause a 
discontinuous behavior of the mass-flux calculations that sometimes led to erratic behavior of the 
convection scheme. This behavior occurred even for vertical configurations that are considered 
to be very good for 3-D models (e.g., layer thicknesses of only ~40 m in the boundary layer and 
the cloud layer). Further work led to changes in the formulations of the shallow convection 
scheme so that the dependency on discrete layer thickness was reduced or eliminated. The 
consequence was that the key quantities of the shallow-convection scheme, including mass flux 
and detrainment rates, now evolve in a generally smooth manner as the clouds develop and grow. 
Not only does the reformulation provide a more realistic shallow-cloud structure, but it also 
makes the evolution of calculated water-budget terms smoother and much easier to interpret. 

Following the extensive testing of each individual parameterization, several combined 
experiments were conducted in the 3-D MM5 using the PLACE or Blackadar Slab land-surface 
scheme, and the improved versions of the PSU 1.5-order TKE-predicting turbulence scheme and 
shallow-convection scheme. These experiments showed that the MM5 can produce fairly 
realistic distributions of shallow convective cloud fields in stratocumulus and fair-weather 
cumulus (humilis) conditions. Cloud area and liquid water content are reasonable and are within 
the accuracy expected, given typical uncertainties in the initial atmospheric and soil states. The 
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vertical structure for the resolved-scale temperature and moisture fields in the cloud layer appear 
to evolve in a generally realistic way, although further testing is required. 

An offline PLACE system forced by observations was constructed and refined for this 
study. This stand-alone version of PLACE appears to produce realistic soil moisture fields 
useful in initializing the 3-D MM5-PLACE model. We have compared the forcings and the 
results of the stand-alone PLACE with available data and obtained reasonable results. We also 
have carried out a test of improved physical parameterizations in the framework of the MM5v3 
mesoscale modeling system. The case chosen was characterized by rather weak synoptic forcing 
typical of summer conditions over the SGP. This July case presented complex mesoscale 
variation in convective and surface conditions and therefore posed a considerable challenge to 
getting realistic modeling results. 

Model results obtained with the new shallow convection parameterization scheme in this 
advanced MM5 system produced realistic spatial patterns of clouds, although cloud amounts 
appear to be low when compared to satellite derived data. We also examined the effects of 
parameterized clouds on the vertical structure of the BL and on the radiation forcing at the 
surface. The preliminary results show physically plausible effects, although more evaluation is 
needed to better ascertain the quantitative skills of the parameterization. Model results that used 
improved initialization and modeling approaches for the surface processes, especially soil 
moisture and vegetation cover, showed that the spatial and temporal structure of moisture in the 
boundary layer is positively affected by the inclusion of the improved surface treatments. In 
particular, the model was able develop an east-west gradient in moisture over central Oklahoma 
that was absent in an experiment using the more standard Slab land-surface scheme. The 
modeled results, however, suffered from sub-grid convective precipitation that was predicted 


198 


over Oklahoma, but did not occur. This precipitation, which was a result of the atmospheric 
initialization, affected in turn evaporation rates in the runs that included the dynamic soil 
moisture parameterization. When this feedback was eliminated (in effect, eliminating the impact 
of the faulty initial thermal stability and moisture fields), the soil moisture values produced by 
the model appeared closer to the observations and the spatial and temporal variation of moisture 
compares better with observations. Moreover, we have completed a preliminary study of the 
effects of considering sub-grid moisture variability in the MM5-PLACE model results. The 
results of these experiments show that it appears to be necessary to include the sub-grid 
variability of soil moisture in order to avoid un-realistic diurnal cycles in evaporation rates, at 
least for the 36-km grid size that we have used. 

More analysis of model results is needed to better evaluate these complex 
parameterization schemes and their interactions. Of special importance appears to be the 
interaction of sub-grid cloudiness, precipitation, radiation and evaporation parameterizations. 
Without adequate interaction of these processes in the model, the results can rapidly diverge 
from the observations. Another point of interest is to determine if the development of a more 
intense shallow cloud field over Oklahoma in nature had an impact in suppressing the convective 
precipitation that developed so rapidly in the model solutions. 

Two major accomplishments were achieved on the basis of the research conducted 
through this NASA fellowship grant. First, a set of advanced parameterizations designed for the 
simulation of the physical processes critical to the atmospheric water budget have been 
developed, tested and refined. These parameterizations for land/atmosphere exchanges, 
turbulence, and shallow convection work together smoothly as a unit. They have been evaluated 
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against observations in a variety of situations and show good agreement with measurements and 
physical expectations. The second accomplishment is that this set of physics, applied in the 3-D 
MM5, was used to show the impact of shallow convection on the vertical structure of 
temperature and mixing ratio over the Great Plains. Future work should apply this advanced 
physics suite to additional cases on the mesoscale, and to different model spatial and temporal 
resolutions and scales. 
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APPENDIX A 


Documentation of Surface Characterization Fields Used in PLACE. (RUN14) 


1 . General remarks 

In general, the MM5V3 parameters have been used when possible. More specialized PLACE 
parameters - mainly related to heterogeneous conditions - have been set at the most simple 
values (non-heterorogeneous). 


2. Soil Parameters 

2.1 BPOWERC 1 

description: 

units for PLACE: 
definition method: 


2.2 THETAS ATC1 

description: 
units for PLACE: 
definition method: 


2.3 PSISATC1 

description: 
units for PLACE: 
definition method: 


2.4 XKSATC1 

description: 
units for PLACE: 
definition method: 


B parameter. It is the exponent in the relationship between 
moisture potential (Psi) and soil moisture (Clapp and Homberger). 
no units 

defined using the table in MM5V3: prmsoi.F and the soil index 
map of MM5V3 


Saturated volumetric soil moisture. 
m3/m3 

defined using the table in MM5V3: prmsoi.F and the soil index 
map of MM5V3 


Saturated soil moisture potential, 
cm 

defined using the table in MM5V3: prmsoi.F and the soil index 
map of MM5V3 


Saturated soil moisture conductivity, 
cm/s 

defined using the table in MM5V3: prmsoi.F and the soil index 
map of MM5V3 
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2.5 Other PLACE soil parameters 


SIGMATHETA 

0.0 

standard deviation in soil moisture content 


COEFFVAR 

0.0 

coeff. of variation for soils (see E&C88 fig. 4) 


TOPFLAG 

0.0 

is a flag used for computation of soil heat flux: 'TOPFLAG' is 1 
when the uppermost layer of the soil is assumed to have a reduced 
conductivity due to loose particles, litter, etc. acting as insulation; 
also due to conditions applying to rapid drying of soil otherwise 
TOPFLAG=1.0. 

— 

DDELTOPO 

0.0 

vertical dimension to calculate horizontal discharge 


DELHORIZ 

50000. horizontal dimension to calculate horizontal discharge 


FRACTURES 1 

1.0 

distribution of cracks and fisures at lower boundary of grid 


MAXSATL1 

5.0 

maximum saturation level: this is a parameter which controls the 
level at which the supersaturation prevention code is activated. 
The prevention of Darcy motion is used from level MAXSATL 
down to the model base for saturated conditions. 

— 

IDSCHG1 

1.0 

horizontal discharge switch for layer 1 


IDSCHG2 

1.0 

horizontal discharge switch for layer 2 


EDSCHG3 

1.0 

horizontal discharge switch for layer 3 


IDSCHG4 

1.0 

horizontal discharge switch for layer 4 


IDSCHG5 

1.0 

horizontal discharge switch for layer 5 

— 

ICOSBY 

0.0 

use Cosby relationships for soil parameters 


NCOLACTV 

1.0 

number of active soil columns (use 1 when ICOSBY=0) 


SCOLWGT1 

1.0 

column weights 

— 

FRSILTC1 

0.0 

silt fraction in columns 


FRSANDC1 

0.0 

sand fractions in columns 


FRCLAYC1 

0.0 

clay fraction in columns 

— 
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3. Vegetation/Land Use parameters 
3.1 ALBEDO 


description: 

Albedo 

units for PLACE: 

no units 

definition method: 

defined using the table in MM5V3: prmveg_usgs.F and the USGS 
map in MM5V3. 

3.2 EMIS 

description: 

Emissivity 

units for PLACE: 

no units 

definition method: 

set to 1 .00 

3.3 ZOHOLD 

description: 

Surface roughness 

units for PLACE: 

cm 

definition method: 

defined using the table in MM5V3: prmveg_usgs.F and the USGS 
map in MM5V3 

3.4 LAI 

description: 

Leaf area index 

units for PLACE: 

no units 

definition method: 

set to 7.00 

3.5 PND VI 

description: 

vegetation fraction 

units for PLACE: 

no units 

definition method: 

defined using monthly climatological maps in MM5V3 

3.6 BM 

description: 

biomass on the surface 

units for PLACE: 

grams of water equivalent / cm2 

definition method: 

defined using educated guess by P. Wetzel and USGS map in 
MM5V3 
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3.7 DEWMAX 


description: 

units for PLACE: 
definition method: 

3.8 PSICRIT 

description: 
units for PLACE: 
definition method: 

3.9 RSUBP 

description: 
units for PLACE: 
definition method: 

3.10 SIGLOGRP 

description: 
units for PLACE: 
definition method: 

3.11 GEMVEG 

description: 

units for PLACE: 
definition method: 


3.12 RSMIN 

description: 
units for PLACE: 
definition method: 


3.13 RSMAX 

description: 
units for PLACE: 
definition method: 


the maximum amount of dew and/or intercepted precipitation on 
the surface 
cm 

set to 0.05 


critical water potential at which stomata close due to water stress 
cm 

set to -2.0E+4 


internal plant resistance between root and leaf 
s 

set to 8.0E+8 


standard deviation of Log of RSUBP/1.0E8 
no units 
set to 0 


species dependent empirical parameter related to the effects of the 
vapor pressure deficit of the atmosphere 

(hPa)' f 


defined using the table in MM5V3: prmveg_usgs.F and the USGS 
map in MM5V3 


minimum stomatal resistance 
same as RSMAX 

defined using the table in MM5V3: prmveg_usgs.F and the USGS 
map in MM5V3 


maximum stomatal resistance 
same as RSMIN 
set to 50. 
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3.14 TEMPSTMN 


description: 
units for PLACE: 
definition method: 

3.15 TEMPSTMX 

description: 
units for PLACE: 
definition method: 

3.16 TEMPSTOP 

description: 
units for PLACE: 
definition method: 

3.17 RGLFLX 

description: 
units for PLACE: 
definition method: 


minimum temperature for stomatal functioning 
K 

set to 273 


maximum temperature for stomatal functioning 
K 

set to 323 


optimum temperature for stomatal functioning 
K 

set to 298 


normalization factor for incoming solar radiation. 

W/m2 

defined using the table in MM5V3: prmveg_usgs.F and the USGS 
map in MM5V3 


3.18 Root fractions at each soil layer 


ROOTFRAC1 0.00 

ROOTFRAC2 0.50 

ROOTFRAC3 0.25 

ROOTFRAC4 0.25 

ROOTFRAC5 0.00 


3.19 Proportionality constant Alphapl 


description: 


proportionality constant between soil-to-root resistance and inverse 
soil hydraulic conductivity 


ALPHAPL 1 0.00 

ALPHAPL2 1.00 

ALPHAPL3 0.03 

ALPHA PL4 0.03 

ALPHAPL5 0.00 
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APPENDIX B 


Observed and Offline PLACE Soil Moisture for the 35 Mesonet Stations with Data 


Note: Each row in the following figures contains the measured and modeled soil moisture time 
series at each Mesonet station with full soil moisture data, at the four levels with measurements. 
Offline PLACE model results have been interpolated to the corresponding measurement levels (5 
cm, 25 cm, 60 cm and 75 cm depth). Ordinates are volumetric soil moisture contents. 

Abscissas are days from July 1997. Lines with crosses correspond to model results, and lines 
with dots are observations. The station number in each figure title is the official number of the 
Mesonet station (see Appendix C for station identification details). 
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0 10 20 30 
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APPENDIX C 


Station List for Oklahoma Mesonet 


STN STN . 
NO. ID 


STATION 

NAME 


LATITUDE LONGITUDE 

Degrees Minutes Degrees Minutes 


Elevation 
(meters ) 


1 

AD AX 

ADA 

34 

2 

ALTU 

ALTUS 

34 

3 

ALVA 

ALVA 

36 

4 

ANTL 

ANTLERS 

34 

5 

ARDM 

ARDMORE 

34 

6 

ARNE 

ARNETT 

36 

7 

BBOW 

BROKEN BOW 

34 

8 

BEAV 

BEAVER 

36 

9 

BESS 

BESSIE 

35 

10 

BIXB 

BIXBY 

35 

11 

BLAC 

BLACKWELL 

36 

12 

BOIS 

BOISE CITY 

36 

13 

BOWL 

BOWLEGS 

35 

14 

BREC 

BRECKENRIDGE 

36 

15 

BRIS 

BRISTOW 

35 

16 

BUFF 

BUFFALO 

36 

17 

BURB 

BURBANK 

36 

18 

BURN 

BURNEYVILLE 

33 

19 

BUTL 

BUTLER 

35 

20 

BYAR 

BYARS 

34 

21 

CALV 

CALVIN 

35 

22 

CAMA 

CAMARGO 

36 

23 

CENT 

CENTRAHOMA 

34 

24 

CHAN 

CHANDLER 

35 

25 

CHER 

CHEROKEE 

36 

26 

CHEY 

CHEYENNE 

35 

27 

CHIC 

CHICKASHA 

35 

28 

CLAR 

CLAREMORE 

36 

29 

CLAY 

CLAYTON 

34 

30 

CLOU 

CLOUDY 

34 

31 

COOK 

COOKSON 

35 

32 

COPA 

COPAN 

36 

33 

DURA 

DURANT 

33 

34 

ELRE 

EL RENO 

35 

35 

ERIC 

ERICK 

35 

36 

EUFA 

EUFAULA 

35 


48 

96 

40 

296 

35 

99 

20 

417 

47 

98 

40 

450 

13 

95 

42 

181 

12 

97 

05 

268 

04 

99 

54 

719 

01 

94 

37 

113 

48 

100 

32 

754 

24 

99 

04 

511 

58 

95 

52 

184 

45 

97 

15 

303 

42 

102 

30 

1267 

10 

96 

38 

281 

25 

97 

42 

352 

47 

96 

21 

241 

50 

99 

38 

559 

38 

96 

49 

301 

54 

97 

16 

228 

35 

99 

16 

513 

51 

97 

00 

347 

00 

96 

20 

237 

02 

99 

21 

587 

37 

96 

20 

209 

39 

96 

48 

291 

45 

98 

22 

361 

33 

99 

44 

692 

02 

97 

55 

329 

19 

95 

38 

213 

39 

95 

20 

186 

13 

95 

15 

221 

41 

94 

51 

303 

55 

95 

53 

248 

55 

96 

19 

200 

33 

98 

02 

421 

12 

99 

48 

604 

18 

95 

39 

203 


220 


LATITUDE LONGITUDE 


STN 

NO. 

STN. 

ID 

STATION 

NAME 

Degrees 

Minutes 

Degrees 

Minutes 

Elevation 
(meters ) 

37 

FAIR 

FAIRVIEW 

36 

16 

98 

30 

407 

38 

FORA 

FORAKER 

36 

50 

96 

26 

330 

39 

FREE 

FREEDOM 

36 

44 

99 

09 

530 

40 

FTCB 

FORT COBB 

35 

09 

98 

28 

420 

41 

GOOD 

GOODWELL 

36 

36 

101 

36 

996 

42 

GRAN 

GRANDFIELD 

34 

14 

98 

44 

342 

43 

GUTH 

GUTHRIE 

35 

51 

97 

29 

328 

44 

HASK 

HASKELL 

35 

45 

95 

38 

183 

45 

HINT 

HINTON 

35 

29 

98 

29 

489 

46 

HOBA 

HOBART 

34 

59 

99 

03 

474 

47 

HOLL 

HOLLIS 

34 

41 

99 

50 

498 

48 

HOOK 

HOOKER 

36 

51 

101 

14 

912 

49 

HUGO 

HUGO 

34 

02 

95 

32 

175 

50 

I DAB 

IDABEL 

33 

50 

94 

53 

110 

51 

JAYX 

JAY 

36 

29 

94 

47 

308 

52 

KENT 

KENTON 

36 

50 

102 

53 

1322 

53 

KETC 

KETCHUM RANCH 

34 

32 

97 

46 

341 

54 

KING 

KINGFISHER 

35 

53 

97 

55 

319 

55 

LAHO 

LAHOMA 

36 

23 

98 

07 

395 

56 

LANE 

LANE 

34 

19 

96 

00 

182 

57 

MADI 

MADILL 

34 

02 

96 

57 

232 

58 

MANG 

MANGUM 

34 

50 

99 

25 

461 

59 

MARE 

MARENA 

36 

04 

97 

13 

331 

60 

MARS 

MARSHALL 

36 

07 

97 

36 

315 

61 

MAYR 

MAY RANCH 

36 

59 

99 

01 

552 

62 

MCAL 

MCALESTER 

34 

53 

95 

47 

230 

63 

MEDF 

MEDFORD 

36 

48 

97 

45 

330 

64 

MEDI 

MEDICINE PARK 

34 

44 

98 

34 

488 

65 

MIAM 

MIAMI 

36 

53 

94 

51 

248 

66 

MINC 

MINCO 

35 

16 

97 

57 

431 

67 

MTHE 

MT HERMAN 

34 

19 

94 

49 

286 

68 

NEWK 

NEWKIRK 

36 

53 

96 

55 

369 

69 

NORM 

NORMAN 

35 

15 

97 

29 

361 

70 

NOWA 

NOWATA 

36 

45 

95 

36 

206 

71 

OILT 

OILTON 

36 

02 

96 

30 

257 

72 

OKEM 

OKEMAH 

35 

26 

96 

16 

263 

73 

OKMU 

OKMULGEE 

35 

35 

95 

55 

205 


221 



STN STN . 
NO. ID 


STATION 

NAME 


LATITUDE 


LONGITUDE 


Degrees Minutes Degrees Minutes 


Elevation 
(meters ) 


PAUL PAULS VALLEY 
PAWN PAWNEE 
PERK PERKINS 
PRYO PRYOR 
PUTN PUTNAM 
REDR RED ROCK 
RETR RETROP 
RING RINGLING 
SALL SALLISAW 
SEIL SEILING 
SHAW SHAWNEE 
SKIA SKIATOOK 
SLAP SLAPOUT 
SPEN SPENCER 
STIG STIGLER 
STIL STILLWATER 
STUA STUART 
SULP SULPHUR 
TAHL TAHLEQUAH 
TALI TALI HINA 
TIPT TIPTON 
TISH TISHOMINGO 
TULL TULLAHASSEE 
VINI VINITA 
WALT WALTERS 
WASH WASHINGTON 
WATO WATONGA 
WAUR WAURIKA 
WEAT WEATHERFORD 
WEBB WEBBERS FALLS 
WEST WESTVILLE 
WILB WILBURTON 
WIST WISTER 
WOOD WOODWARD 
WYNO WYNONA 
NINN NINNEKAH 
ACME ACME 
APAC APACHE 


34 

43 

97 

14 

292 

36 

22 

96 

46 

283 

36 

00 

97 

03 

292 

36 

22 

95 

16 

201 

35 

54 

98 

58 

589 

36 

21 

97 

09 

293 

35 

07 

99 

22 

538 

34 

12 

97 

35 

280 

35 

26 

94 

48 

156 

36 

11 

99 

02 

540 

35 

22 

96 

57 

329 

36 

25 

96 

02 

285 

36 

36 

100 

16 

771 

35 

33 

97 

20 

373 

35 

16 

95 

11 

175 

36 

07 

97 

06 

272 

34 

53 

96 

04 

254 

34 

34 

96 

57 

320 

35 

58 

94 

59 

290 

34 

43 

95 

01 

204 

34 

26 

99 

08 

388 

34 

20 

96 

41 

268 

35 

50 

95 

25 

189 

36 

46 

95 

13 

236 

34 

22 

98 

19 

308 

34 

59 

97 

31 

340 

35 

51 

98 

32 

516 

34 

10 

97 

59 

285 

35 

30 

98 

47 

538 

35 

28 

95 

08 

146 

36 

01 

94 

39 

348 

34 

54 

95 

21 

201 

34 

59 

94 

41 

143 

36 

25 

99 

25 

625 

36 

31 

96 

21 

269 

34 

58 

97 

57 

358 

34 

48 

98 

00 

397 

34 

55 

98 

17 

440 


